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Grant NCR 43-003-015. It is the masters research of the first 
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INTRODUCTION AND LITERATURE REVIEW 

During steady operation of a liquid propellant rocket engine the 
injected propellants are converted hy various physical and chemical 
processes into hot burned gases which are subsequently accelerated to 
supersonic velocity by passing through a c on verg ing -d i verg ing nozzle. The 
operation of such an engine* however* is seldom perfectly smooth. Instead 
the quantities which describe the conditions inside the combustor (i.e. 
pressure, density, temperature, etc.) are time-dependent and oscillatory. 
Such oscillations can be of either a destructive or nondestructive nature. 
Nondestructive unsteadiness is characterized by random fluctuations in the 
flow properties and includes the phenomena of turbulence and combustion 
noise. Unsteady operation of a destructive nature, on the other hand, is 
characterized by organized oscillations in which there is a definite 
correlation between the fluctuations at two different locations in the 
combustor. Such oscillations have a definite frequency and result in 
additional thermal and mechanical loads that the system must withstand. 

Unsteady operation of the destructive variety, knovn as combustion 
instability, was first encountered in 1940. At that time a British group 
testing a small solid-propellant rocket motor observed sudden increases 
of pressure to twice the expected level, enough to destroy a motor of 
flight weight. Since tret time every major rocket development program 
has been plagued by combustion instability of some form. These 
oscillations in the combustion chamber con have several detrimental effects. 
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la soTOt casts, particularly in solid-propellant rockets , testability 
can causa the steady-stats pressure to increase to a point at which the 
rocket motor will explode. In liquid-propellant rocket chambers experi- 
encing unstable combustion, heat transfer rates to the walls considerably 
exceed the corresponding steady state heat transfer rates, resulting in 
burn-out of the walls. If the chamber can survive these effects, mechanical 
vibrations In the rocket system can cause mechanical failure or destroy the 
effectiveness of the delicate control and guidance systems. 

The phenomenon of combustion testability depends heavily upon the 
unsteady behavior of the combustion process. The organized oscillations of 
the gas within the chamber must be coupled with the combustion process in 
such a way as to form a feedback loop. In this manner part of the energy 
stored in the propellants becomes available to drive large amplitude 
oscillations. An understanding of this coupling between the combustion 
process and the wave motion is necessary in order to predict the stability 
characteristics of rocket engines. 

Combustion instability problems in liquid propellant rocket motors 
usually fall into one of three categories according to the frequency of 
oscillation. Low frequency combustion instability, also known as chugging, 
is characterized by frequencies ranging from ten to several hundred 
hertz, nearly spatially, uniform properties, and coupling with the feed 
system of the rocket. This type of instability is less detrimental than 
other forms, and the means of preventing it are well understood. Low 
frequency instability will not be considered. 

A second type of combustion instability, which is less frequently 
observed, has a frequency of several hundred cycles per second. This 
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type of oscillation la associated with the appearance of*entropy waves 
inside the combustion chamber. 

The third and aost important form of combustion instability is 
known as high frequency or acoustic instability. As the nano suggests , 
this type of instability represents the case of forced oscillations of the 
combustion chamber gases which art driven by the unsteady combusticn process 
and interact with the resonance properties of the combustor geometry. The 
observed frequencies* which are as high as 10*000 cycles per second* are 
very close to those of the natural acoustic modes of a dosed-ended 
chamber of the same geometry as the me experiencing unstable combustion. 
High frequency combustion instability is by far the most destructive and 
is the type to be considered by the following analysis. 

High frequency combustion instability can resemble any of the 
following acoustic modes: (1) longitudinal, (2) transverse* and (3) 

combined longitudinal-transverse modes. Longitudinal oscillations are 
usually observed in chambers whose length to diameter ratio is much greater 
than one; in this case the velocity fluctuations are parallel to the axis 
of the chamber and the disturbances depend only on one space dimension. 

For much shorter chambers the transverse mode of instability is most 
frequently observed. Transverse oscillations in rocket motors are 
characterized by a component of the velocity-perturbation which is 
perpendicular to the axis of the chamber but the disturbances can depend 
upon three space dimensions. Such oscillations can take either of two 
forms: (1) the standing form in which the nodal surfaces are stationary 

and (2) the spinning form in which the nodal surfaces rotate in either the 
clockwise or counterclockwise direction. Transverse combustion insta- 
bility, particularly that resembling the first tangential mode* has been 


Since the early 1950's much experimental and analytical research 
has been devoted to better understanding the phenomenon of high frequency 
combustion instability. Most of the theories presented prior to 1966 were 
restricted to circumstances in which the amplitudes of the pressure 
oscillations were infinitesimally small in the linear regime. Prominent 
among these are the pioneering studies of longitudinal instability by 
Croc co [1] as well as the studies of transverse instability by Scala [2], 
Reardon [3], and Culick [43. A complete discussion of these theories, is 
given in the work of Zinn [5] and will not be repeated here. 

Although linear theories provide the propulsion engineer with 
considerable insight into the problem, their applicability and usefulness 
in design is limited. The linear theories cannot provide answers to such 
important problems as the limiting value of the pressure amplitude 
attained by a small disturbance in the case of a linearly unstable engine, 
or the effect of a finite-amplitude disturbance upon the behavior of a 
linearly stable engine. In the latter case the result of many tests 
indicate that undfer certain conditions the introduction of sufficiently 
large disturbances into a linearly stable engine can trigger combustion 
instability. Another shortcoming of linear theories is the fact that 
their predictions cannot be compared directly with available experimental 
data; for, in the majority of cases, the experimental data is obtained 
under conditions in which the combustion instability is fully developed 
and in a non-linear regime. Therefore, theories accounting for these 



nonlinearities associated with combustion instability are needed. A 
more detailed discussion of the nonlinear aspects of combustion instability 
can be found in a work by Zinn [5]. 

In the field of finite amplitude (nonlinear) combustion instability, 
mathematical difficulties have precluded any exact solutions, and 
approximate methods and numerical analysis have been used almost exclusively. 
For this reason publications in this field are scarce. Notable among these 
is the work of Maslen and Moore [6] who studied the behavior of finite 
amplitude transverse waves in a circular cylinder. Their major conclusion 
was that, unlike longitudinal oscillations, transverse waves do not steepen 
to form shock waves. Maslen and Moore, however, considered only fluid 
mechanical effects; they did not consider the influences of the combustion 
process, the steady state flow, and the nozzle which are so important 'in 
the analysis of combustion instability problems. Nevertheless, pressure 
recordings taken from engines experiencing transverse instability reveal 
the presence of continuous pressure waves similar in form to those 
predicted by Maslen and Moore. 

One of the first nonlinear analyses to include the effects of 
the combustion process and the resulting steady state flow was performed 
by Priem and Guentert £7], In this investigation, the problem was made 
one-dimensional .by considering the behavior of tangential waves traveling 

t 

in a narrow annular combustor of a liquid propellant rocket motor. They 
used a computer to solve numerically the resulting nonlinear equations for 
various values of the parameters involved. Due to the many assumptions 
involved in the derivation of the one-dimensional equations, the results 
of this investigation are open to question. 
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The successful use of the time-lag concept (see Crocco [1]) in the 
linear theories prompted a number of researchers to apply this model to 
the analysis of non-linear combustion instability. By considering a 
chamber with a concentrated combustion zone and a short nozzle, Sirignano 
[8] demonstrated the existance of continuous, finite-amplitude, longitudinal 
periodic waves. These solutions were shown to be unstable, however , thus 
indicating the possibility of triggering longitudinal oscillations. 

Mitchell (9] extended the work of Sirignano to include the possibility of 
discontinuous solutions. In this manner he was able to show that the final 
form of triggered longitudinal instability consisted of shock waves moving 
back and forth along the combustion chamber. Mitchell also considered the 
more realistic case of distributed combustion. 

In the analyses of Priem, Sirignano, and Mitchell the oscillations 
were dependent on only one space dimension. One of the first researchers 
to study finite-amplitude three-dimensional combustion oscillations was 
Zinn [5] whose work is an extension of the linear transverse theories and 
the analysis of Maslen and Moore. Using Crocco* s time lag model Zinn 
investigated the nonlinear behavior of transverse ' waves' in a chamber with 
a concentrated combustion zone at the injector end and an arbitrary 
converging-diverging nozzle at the other end. In this case, it was 
necessary to ext-end Crocco* s burning rate expression and transverse nozzle 
admittance relation to obtain the appropriate boundary conditions for the 
case when the flow oscillations are of finite size. As a result of this 
analysis Zinn was able to prove the existance of three dimensional 
finite-amplitude continuous waves which are periodic in time. In 
addition, he was able to prove the possibility of triggering combustion 
oscillations. An analytical criterion for the determination of the 
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stability of such waves was derived, but because of its complicated form 
and the limited capacity of available computers no specific numerical 
results were obtained. 

In more recent years other investigators such as Burstein [10] 
have attempted to solve numerically the equations describing instabilities 
that depend on two space dimensions. Although the resulting solutions 
resemble experimentally observed combustion instability, this method 
requires excessive computer time, and studies of this type for three- 
dimensional oscillations will have to await the development of a much 
faster breed of computers. 

In a recent publication by Powell [11], the problem of analytically \ 
and numerically analyzing multidimensional non-linear combustion instability 
was investigated. The problem in doing this is that a system of non- 
linear coupled partial differential equations whose solutions must 
satisfy a complicated set of boundary conditions governs the phenomena of 
combustion instability. These boundary conditions may describe the 
unsteady burning process of the wall of a solid propellant rocket motor; 
the conditions at an idealized concentrated combustion zone of a liquid- 
propellant rocket engine; or the unsteady flow of the entrance of a 
converging-diverging nozzle. Previously, in an effort to obtain analytical 
solutions to various combustion instability problems, investigators have 
been forced to simplify the original problem to such an extent that it no 
longer resembled the real problem that originally was to be solved. Powell 
proposed a method to perform a nonlinear stability analysis with relative 
ease. This methodj applicable to both linear and non linear problems with 
complicated boundary conditions, was a modified form of the classical 
Galerkin method. The Galerkin method [11] is an approximate mathematical 




technique vhich has been successfully employed in the solution of various 
engineering problems in the field of acoustics. Powell used this method 
to specifically study the non-linear behavior of combustion driven 
oscillations in cylindrical combustion chambers in which the liquid 
propellants are injected uniformly across the injector face and the 
combustion process is distributed throughout the combustion chamber. Based 
upon the results of his second and third order theories, the following 
nonlinear mechanisms were found to be important in determining the non- 
linear stability characteristics of the system: (1) the transfer of energy 

between modes, (2) the self -coupling of a mode with itself and (3) a non- 
linear combustion mass source. Powell found that the self -coupling 
mechanism was important in the initiation of triggered instability, while 
the non-linear driving mechanism was important in the determination of the 
final amplitude of triggered instability. 


Statement of the Problem 


In this thesis, the problem of velocity-sensitive instability will 


be considered. Based upon previous work on this problem, only transverse 
oscillations will be considered due to mathematical simplicities. Also, 
the specific geometry of the combustion chamber to be analyzed will be 
annular or ring-like. The purpose of this thesis is to investigate the 
mechanisms which cause these instabilities due to the combustion process 
in a liquid propellant annular combustion chamber and attempt to state 
which mechanisms or conditions impose the greatest effect upon stability 
of combustion. 

In Chapter 2 of this thesis, the governing equations of fluid 
motion (i.e., balance of mass and momentum) are stated. From the equations 
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the general acoustic wave equation for non-linear combustion is derived. 

In this derivation, both steady state and deviations from the steady-state 
conditions are considered and their effects incorporated into the general 
acoustic wave equation. 

In Chapter 3, the Galerkin method is used to obtain, from the 
general acoustic equation of Chapter 2, equations governing the modal 
amplitudes associated with the first two modes of transverse oscillation 
in a thin annular combustion chamber. These equations for the annular 
combustion chamber are solved numerically by the use of a Runge-Kutta 
program for various conditions. 

In Chapter 4, a set of approximate equations are derived from the 
modal amplitude equations presented in Chapter 3 by use of the two-variable 
perturbation technique. These resulting approximate equations are 
expressed both in the modal amplitude and amplitude-phase angle form. In 
this chapter, four special cases are presented for which closed-form 
solutions can be found. These four cases are (1) standing wave — no 
combustion, (2) standing wave — no gas dynamic nonlinearities, (3) 
traveling wave — no combustion, and (4) traveling wave — no gas dynamic 
nonlinearities. For problems not falling within the above categories, 
a numerical analysis is employed to solve approximate equations. 

In Chapter 5, the results contained in the previous two chapters 
are discussed and compared. Stability limits are obtained and the effect 
of neglecting various physical effects are discussed. In addition, the 
accuracy of the perturbation method is evaluated. A summary o i the 
research contained in this thesis is presented in this chapter. 

In Chapter 6, a statement of conclusions is made along with 
recommendations for future research in this area. 


Chapter 2 


DERIVATION OF THE GOVERNING ACOUSTIC HAVE EQUATION 

In order to investigate the non-linear combustion instabilities 
that occur in liquid propellant rocket engines , one must start with the 
balance laws of mass and momentum. Also* for this problem, a constitutive 
equation was formulated relating pressure and density. Mathematically, 
these principles are respectively 


3p* 

3t* 


ft 
-*> 
+ V 


ft * 

(p u) = B* 


( 2 . 1 ) 


Hat* + u 



( 2 . 2 ) 


ft *2 * 
p = a p. 


(2.3) 


where 

ft 

p - gas density 

ft 

t - time 
* 

t 

ft . 

V - del operator of the system i + j + K 

ft 

u - velocity of the gas 

it 

B - fuel drop burning rate per unit volume 
ft 

p - pressure of the gas 
*2 

a - constant of proportionality (in this case - speed of 
sound ) . 
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The * representation denotes that the above physical quantities are 
dimensional. Equations (2.1) - (2.3) are based on the assumption that 
the fuel drops serve only as a source of mass for the gas phase. 
Interphase transfer of momentum and energy are neglected. 

Combining equations (2.2) and (2.3) , the resulting equation is 


u ■* * j *2i* 

'V 3 t* + u * * uj = - a Op. 


For the physical situation depicted in Figure 1 


inlet 


gas- liquid drop 
mixture _ 


exhaust 


combustion chamber 


variable area cross section 


fuel drops enter here 
through injector plates 


Figure 1. Schematic of a Liquid Propellant Combustion Chamber 
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A convenient non-dimensionalization of the variables is as follows: 


e e ,e 

P * P Q P (P Q - initial density of gas) 


* * 

■+ * + 

u = a u 


$ * 4-$ 


L 


. L* 

** 8 i?T t 


*2 * 

P* * a p Q p 


ft ft 

5 - ° 

B “ 1JT 


B. 


Substituting these non-dimensional relations into equations (2.1), (2.3), 
and (2.4), the results are 


!r + $ • (pu) * B 


(2.5) 


P 


3u 

at 


+ p u 


$ u 



( 2 . 6 ) 


p = P 


(2.7) 


where the unstarred quantities are dimensionless. 




Dividing through by density p, equation (2.6) becomes 


3u * ie vo 

u * ?u * ’ p 


Since , 


I s - 8 $ ** P. 


the governing equations can be summarized as 


|| + t • (pu) * B 


l^+u • ? u * -ttn.fi 


P = p- 


It will now be shown that to the order of approximation inherent 
in these equations, the flow is irrotational , that -is t. x u * 0. To do 
this, take the curl of equation (2.10) and set it equal to zero. The 
resulting equation becomes 


# 


Since the curl of any gradient is zero. This may be rewritten as 


?x |^ + ? x (u * K) ! 0, 


(2.13) 
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; The vopticity ft is defined to be 

2 * ? x u. (2.14) 

Thus, 

^ x !lF *ft ( ^ x “ )a ||’ (2,15) 

From the vector identity 

I • ? X * ?(ljP) - I X <* X X) 

it follows that 

u . $ u « $(* 5 *^) ■ u x2. (2.16) 

Therefore , 

ft x (u • ft u) = ft x [ft(lju 2 ) - u .*.$] (2.17) 

Recognizing that the curl of any gradient is zero, equation (2.17) 
reduces to 

ft x (u • ? I) : • ? x (u x 2) , (2.18) 

Using the vector identity 



7 x (t x h = (ft • ft)ft - ft(ft • ft) - (ft • ft)ft + ft(ft • I) 
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•quation (2.11) can be expressed as 

3 x (u • 3 u) * - t(3 • 3) u - 3(3 • u) -‘Xu • 3)3 

♦ u(3*3)]. (2.19) 

Therefore, equation (2.13) becomes 

|| - (2 • ^)u ♦ 3(3 • u) + (u • 3)3 - u(3 • 3) « 0. (2.20) 

»t 

Equation (2.20) can now oe modified by using the definition for the total 
(comoving) derivative which is 

* • * - • 

Substituting this expression into equation (2.20) and simplifying, the 
resulting equation becomes 

j|. 3 -(3 3> - (3 3). u ♦ 3 (3 • 3). (2.2D 

Rewriting u (3 - t 3) as u [3 • (3 x u)] which is zero since the divergence 

0 

of the curl of any vector is zero, equation (2.21) becomes 

4 

{§. S.(f u) -(«*>• S . (2.22) 

Dt 

The implications of this equation for a fluid starting from rest are as 
follows. At the initial instant of time (t * 0), the vorticity of any 
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fluid particle will be zero. Thus, ths time dsrivativs of ths vorticity 
of ths part ids will bs zero, implying that ^ f 0 « t * 0, Sines 
2 * 0 and * 0 at t ■ 0, it follows that 8 ■ 0 at ths nsxt instant of 
tins. By induction, it can bo shown that 2 * 0 for all tins unlsss ths 
velocity gradient becomes infinite for any t ■ 0. It is sssunsd in what 
follows that this doss not occur wad ths flow is treated as irrotational. 

Since izrotationality has been proven, ths velocity vector u can 
bs expressed as 


U ■ 2 f (2.23) 

where ifr is the velocity potential. Substituting equation (2.23) into the 
left hand side of equation (2.10), the result is . . 

, ; . ? ; . || . . Z * 3 

* §J- (3 *> ♦ 3[>i(? *> J ] - ? * K 3 . 12.2k) 

For irrotational flow ( 2 * 0), the right hand side of equation (2.24) 
becomes 

* 

* 

+ hi* * ^ # C2.25) 

Therefore, equation (2.10) can be written as 


2[Jj- t • $*> ♦ * ■ P « 0 . 


(2.26) 
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Spatially integrating aquation (2.26) produces 


|i t H ♦ In p * e(t) 


where a(t) is a function of integration. From equation (2.23), it can 
be seen that an arbitrary function of time can be added to $ without 
affecting the result for u. Thus, a(t) could be absorbed into ♦. The 
same thing is accomplished by setting a * 0 which results in 


tn p * - - *$♦ • *\> 




Thus, p and u are both known as functions of $>, From equation (2.9), the 
governing equation for $ can be written symbolically as 


|| + p? 2 i|» + % • % * B 


"(ft * ** 

P « P * e ' ' . 


(2.30.b) 


Rather than combining these quantities immediately, it is convenient to 
first make further simplifications based on the nature of the physical 
problem that it is desired to analyze. 
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Steady State Solution 

First, the steady state solution of equations (2.30) corresponding 
to purely axial motion will be found. Define the steady-state velocity 
potential $ by 

<|» = e<j>(z) (2.31) 

where e (assumed small) is the measure of the deviation of the density 
from its initial value (see equation 2.32 below). The bar notation will 
represent steady-state conditions. The steady-state burning rate u is 
defined from 


B = i(z). - (2.-32) 

While many other situations are possible, attention will be confined in 
the present work to the case when u s 0(e). To indicate this let 

53 = eo (o = 0(1)) . (2.33) 

Thus, the burning rate B can be expressed as 


B = eo. 


t 

Equation (2.30.b) can now be written 



(2.34) 


(2.35) 


Using the Taylor series expansion for the exponential function and 
retaining only the first two terms, equation (2.35) becomes 


p = 1 - ^n 2 


+ . . 
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(2.36) 


Substituting equations (2. 31), (2.34), and (2.36) into equation (2. 30. a) 
and dividing the result by e yields 


or 



(2.37) 


d 2 $ 3 

diT " 2 



(2.38) 


Retaining only terms of 0(1) produces 


d 2 $ 

diT 


o. 


(2.39) 


For simplicity, only the case of uniformly distributed combustion (i.e. 
o = constant) will be considered. Thus, integrating' equation (2.39) one 
obtains 


5 

dz 


o z + C. 


where = u is the steady state velocity of the gas. 

At the injector (Z = 0), u = 0. Thus, = 0 and 



(2.40) 


u 


0 z . 


(2.41) 


Deviations from Steady State 

It is now desired to investigate the stability of the steady 
state solution discussed above. Toward this end, an additional 
velocity potential related to perturbations from the steady state is 
defined by the equation 

<1/ = c Ci + $(x, y, z, t )]. (2.42) 

A perturbation burning rate B is also defined by the equation 

B = u> + ew. (2.43) 

It is assumed that u = 0(e) and this is indicated by defining a. function 
o such that a = 0(1) and w = oe. Then equation (2.43) becomes 

B = e(o + eo) . (2.44) 

Taking the gradient of equation (2.42), one obtains 

7\ji = e[V $ + ^ (2.45) 

or 

t 

% = e[u e + 7 <j>] . (2.46) 

z 


From equation (2.42), the time derivative of ifj can be expressed as 


Substituting the equations (2.46) and (2.47) into equation (2.30.b) 
and simplifying, one obtains 



(2.48) 


3 

Expanding (2.48) in a Taylor series and neglecting terms of 0(e ) and 
higher produces the expression 


P = P 


1 


e + e 2 [-*G 2 + $<|> * ^$) 



(2.49) 


Substituting equations (2.42), (2.44), and (2.48) into equation (2. 30. a) 


dividing 

the result by t 

: leads to 



"9^ + 

■ [-* k ( " * 

^ * ^$) 


+ 

. . . + Jl - e 

£ +e2 [ 

-? 5 (u 2 + ^ • 7$) - 

- d<t> 
U 3z 

+ 

•(§* 

^ 2 $j + 

(ue z +7$) • eo 

9$ 

It 

+ 

e 2 7(- i j(u 2 + ^ 
* 

• - u 

& * h {W] = 

5 + ea 


Neglecting all terms of 
state solution that u = 


0(e 2 ) and higher and recalling from the steady- 
— = oz and 4^ = o yields 




It 5 * “ lift + (*♦ ‘ * It)] = *° e * 


Substituting 


% »_ (?*.?♦) = % . ? a. 


(2.52) 


into equation (2.51), results in 


0 - V2* t e [ 2 ($* • ? |f) + 2 3 


. at - 

3z3t 3t 


fc (»*♦ - 0 )] ■ “ 


(2.53) 


where only terms of 0(1) and 0 (e) have been retained. Equation (2.53) can 

02a 

be further simplified by observing that 7 2 <j» = -s-*- + 0 (e). 


Thus, the last term of equation (2.53) can be written 




Since the other terms of 0(e 2 ) have already been neglected, consistency 
requires that this term be deleted and the equation be rewritten as 


- V 2 $ + 


[ 2 ( ? * • 


+ 2u + 


H 5 ]-- 


(2.54) 
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Ht this thesis, attention will be confined to transverse instability. 

For this situation 

* = *(x, y, t). (2.55) 


Therefore, equation (2.54) becomes 

|^ - 72* + ej^V* . $ fi) 

To account approximately for frequency changes due to baffles, nozzle 
shapes, etc., a correction term of the form 





at . 


(2.56) 


was introduced into equation (2.56). This form, one of many possible, was 
chosen so that the linearized form of equation (2.56) would reduce to Love's 
equation for a one-dimensional problem. This linearized form of (2.56) is 

lit _ lit _ £ vJX = o 

3t? 3x2 e K 3x23t2 (2.58) 

Thus, it can be seen that the value K will affect the acoustic frequencies. 

Physically , this is the purpose of baffles , nozzle shapes , and other 

« 

physical parts of the combustion chamber. Therefore, inserting the 
correction term into equation (2.56), the resulting equation becomes 


(2.59) 


where K is ths correction factor. This non-linear wave equation will be 
the basis for numerically and analytically investigating the transverse 
combustion stability problems occurring in liquid propellant rocket 
engines . 


Chapter 3 


DERIVATION OF WAVE EQUATIONS BASED UPON AN 
ANNULAR COMBUSTION CHAMBER 

In Chapter 2, there were no restrictions concerning the geometry 
of the combustion chamber in the derivation of the acoustic wave ] 

equation. In this chapter, however, a set of equations will be developed 
based upon a narrow annular combustion chamber. A typical cross-section 
for such a combustion chamber is shown in Figure 2 below in dimensional 
and dimensionless form. 



(a) Dimensional (b) Dimensionless 

Figure 2.* Dimensional and Dimensionless Form of a Circular 
Cylindrical Combustion Chamber 


* 

In Figure 2 (a), the dimensional quantities are 

ft 

r - radius of a typical point in the combustion chamber 
ft 

R - inside radius of the combustion chamber 
a 

b - thickness of combustion chamber's cross-section. 



In Figure 2 (b), the dimensionless quantities are 

r - non-dimensional radius of a typical point 



1 


The first major assumption to be made in the geometry of the combustion 
chamber is 


|«i (3 

which states that the circular cylinder can be thought of as a thin 
(ring-like) annulus. 

Define the characteristic length L* by 

L* = R*. (3, 


In restricting the analysis to an annulus, a transformation to polar 
coordinates is convenient. Recall that the gradient and Laplacian 
operators in polar coordinates are 


L + t *- 

e r 3r r 36 z 3z 


V 2 


3 2 1 3 . 1 3 2 . 3 2 

" 3r^ + r 3r ^ 36? + 3?T * 


( 3 . 


The second major assumption for the simplification of the velocity 
potential is restricting 
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♦ ■ *<e, t) 


r 4 i . 


(3.4) 


Therefore, using the operators of equations (3.3) on the function of 
equation (3.4) t the results are 


h 


e ae 


= 



(3.5) 


Substituting the results of equation (3.5) into the general acoustic wave 
equation (2.58), the modified wave equation becomes 


3 2 <J> 

It 7 


3 2 $ 

H 7 


+ c 



2 


3 £ 3 2 4 » 

36 * 363t 


K 


is 1 

3t*36*J 


at , 


(3.6) 


Now, express the velocity vector 

u = u + u' (3.7) 

* 

4 

where u - steady-state velocity vector 
u' - perturbation velocity vector . 

From the steady state solution in Chapter 2, the velocity vector was 


defined as 
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Define the perturbation velocity vector by 

& i Si (s>9) 

Substituting equation (3,8) and (3.9) into equation (3.7) and using 
equation (2.23) results in 

« ■ e + e Jf « e ■ (3 - 10) 

To determine only the transverse velocity component of the perturbation 
velocity vector , subtract the perturbed velocity component along the 
axial (z) direction of the chamber from the total perturbation velocity 
vector. Thus, 



■*! 

U 


u e . 
z z 


(3.11) 


In this case, since u = u(6, t) only, there is no perturbed velocity 
component in the axial direction; therefore. 



e 




(3.12) 


It is now desired to find the burning rate 0 in terms of the parameters 
in the wave equation. To obtain this expression, assume velocity sensitive 
combustion with no history effects. Mathematically, the burning-rate 
function for velocity-sensitive combustion will be expressed by the purely 
phenomenological equation 
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Then, the burning rate can be expressed as 



(3.16) 


where j * 0 - no time delay 

* 

1 - time del&y . 

Therefore, substituting equation (3.16) into equation (3.6), equation 
(3.6) can be rewritten 
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(3.17) 


Thare is no closed form solution of equation (3.17) that appears likely. 

The main purpose of the present work is to determine the modifications of 
solutions of the usual acoustic wave equations that are caused by the 
presense of the nonlinear terms multiplied by c in equation (3.17). 

Thus, rsther than attempt a finite di ierence numerical solution of 

i 

equation (3.17), the following procedure was adopted. 

The solution is represented by the Fourier series 

♦(0, t) * fj(t) cos 0 ♦ fj(t) cos 20 + g^(t) sin 0 

+ gj(t) sin 29 ■*• . . . (3,18) 

and initial conditions are chosen such that in the absence of the nonlinear 

terms, the exact solution can be formed using only the first two terms of 

the Fourier series. Because of the quadratic nature of the non-linearities, 

0 

the second two terms in equation (3.18) represent a complete first order 
correction to the acoustic solution due to non-linear gas-dynamic and 
combustion effects. Only the first four terms in equation (3,18) are, 
therefore, retained and the approximate solution determined by this method 
is the simplest one capable of illustrating the influence of the nonlinear 
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terms. The approximation can, of eourss, bs improved by retaining 
additional terms in equation (3.18) but this is not investigated. 

Substituting equation (3.13) into equation (3.17) and using 
the multiple angle formulas t.> simplify terms containing products of 
trignomet&c functions, one obtains * 



sin 20 ♦ . . . ■ 0. 


(3. IS) 
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Equation (3.19) is a summation of terms composed of some function of 
time t and a term containing 0 variation. Since the equation must be 
valid for all values of 0, each of the time dependent coefficients 
of the 9- terms must individually be equal to zero. Therefore, four 
ordinary differential equations governing the time-dependent modal 
amplitudes f^ g^, f 2 and g 2 emmerge from this analysis as the governing 
equations to be used for analysis of instability in an annular combustion 
chamber. Th-..e equations are 


d 2 f 

dtZ 1 + f l 


df 


+ “ dt + 26 


5 

dt 


+ f 


5 

1 dt 


* 

dgj dg 2 

+ g 2 dt“ + g l dt”. 



g l g 2j 


2jenu 


[ 


f It” 2t + g lT g 2T| 


r] 


(3. 20. a) 



dg 


df 

dg 

df 

♦“dt 1 

+ 2e 

g 2 dt 1 

' + f l dt 

g l • dt -• 


dg ] 
f — 1 

2 dt 



f 2 g l 


J - 2jeun £l 


f lT g 2T " f 2T g lT| 


= 0 
(3.20.b) 



+ u 



+ e 


dg 
g l dt 



+ 4Ke 



• *1 ■ 
+ heun g a 2 - f 1 2 - *sjeun g^ - 


0 


(3.20.c) 


u 


d 2 g dg 

dt? + 4 ®2 + 5 d? 


df 

*iai L+f 





t jnwe 



0 . 
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(3.20.d) 


to the following work only instantaneous combustion will be considered. 
Thus, the appropriate equations are equations (3.20)' with j - 0. These 
equations are recapitulated below. 




+ 2c 



+ f 


df 
—2 
1 dt 


+ g 


dg 
2 dt 


♦ g 



1 dt 


+ Ke 


d 2 f 

dt^ 


2neu ^f x f 2 + g x g J * 0 


(3. 21. a) 


d 2 g _ dg 

dt? + g l + “ dt 


+ 2e 


df 


5 2 dt 


+ f 


5a 

1 dt 


df 


n dt 


dg 

- f — J LI 
2 dt 1 



J = 0 (3.21.b) 




(3.21.c) 
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d 2 g„ _ dg, 

3 ^ + “ g 2 ♦ “ at^ 


df 

'1 dt + f 




(3.21-d) 


The equations of (3.21) were solved numerically by the use of the 
quart ic (fourth-order) Runge-Kutta method. To use this method, the 
equations of (3.21) are modified by defining the quantities 


df 



dg ? 

dt = b 2 . (3.22) 

Substituting these expressions into equations (3.20) and solving these 

equations for the highest derivative (in this case - second order), we get 

« 

dt 1 ' [" f l ' “ (a l > " 2Ei ( f 2 (a l } + f l (a 2 } + «2 ( V 

+ gjCbj)) - 2new(f 1 f 2 + g^)! /(I + Ke) 
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_L = [L gi . 5( bl ) - 2ei(g 2 (a 1 ) ♦ f^bj) - g^) 


- f 2 ( bl )j - 2neui(f^g 2 - f^)] 


/(I + Ke) 


dr = [- 4f 2 - “ (a 2> “ ei ( g l ( V 1 f l (a l ) ) 


-Vwnfg^ - f^ 


/(I + 4Ke) 


db 

2 

dt 


[- 


4g. - o)(b ) + ei 


1 ( 8 i (a i > ♦ w) 


+ urns /(l + 41fe) 


(3.23) 


where i is the gas-dynamic index. 

By the development of a computer program incorporating the Runge- 
Kutta algorithm which can solve systems of first-order ordinary differen- 
tial equations, the eight equations (3.22) and (3.23) were numerically 
solved for the eight variables a 1 , a 2 , b^ b 2 , f 1 , f 2 , g^ and g 2 - 
Different cases involving varying the gas-dynamic index, interaction 
index, the correction variable (K), and the order term (epsilon) will be 
discussed and compared with the perturbation method of solution in a 
later chapter. In Appendix B, a sample program listing this calculation 


appears . 


Chapter 4 


TWO-VARIABLE PERTURBATION METHOD APfLIED TO THE 
ACOUSTIC WAVE EQUATIONS 

In this chapter, a set of approximate equations will be developed 
from the governing equations for the modal amplitudes (3.21), by the use 
of the two-variable perturbation method. The two-variable method is well 
suited to this type problem since one expects the solution to consist of 
sinusoidal functions with slowly varying amplitude. Applying this method, 
define two variables representing time 

5 = t 

n = et . (4.1) 

Therefore, the four modal amplitudes would now be 


f ± = fjtt.n) 

f 2 = f 2 (e.n) 

g i = SjU.n) 

g 2 = g 2 U,n) . (**-2) 

By applying the chain rule of differentiation, it can be shown that 


d_Z 3_Z 3Z 

dt " 35 3n 


(4.3) 


and 
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dt 2 H 2 


where Z 3 f^, f 2 , g^, g 2 respectively for each of the above equations. 
By substituting equations (4.3) and (4.4) for each modal amplitude into 
equations (3.21) and keeping terms only of 0(1) and 0(e), the resulting 
equations become 


2 f, 3 2 f, 3f, 3f, 3g, 

TT- + + e[2 — r + + 2f 0 + 2f, -r/ + 2g„ 


+ 2n< J) (f.f 0 + g.g„)] = 0 


+ g l 

+ 6[2 3l3^ + + 

2g 2 + 

2f l 

3 g l 

3 2g l 



n 

+ K -j£- + 2n (1) (f 1 g 2 

“ f 2 g l^ 

= 0 

a. ti 

32f 3f 

_ 3S 1 

_ £ 

+ 4f 2 

3C3n 35 

g l 35 


3 2f 2 




5F 1 

+ *£ ~n(g 1 2-f 1 2)] = 

0 



3 2 g o _ 3g 0 

3f x 


+ 4g 2 

+ el - 2 3C3n + 0 IT - 

“ g l 35 

' f l 

g 3 

7T~ " 

wn(f.g,)] = 0 . 







From the straight-forward perturbation method, define the modal amplitudes 
by the series expansions 

* f 10 (£»h) + c f^e.n) + . . . 

f 2 a f 2o^* n ^ + e f 2i^ ,n ^ ■*■••• 

g l * 8 10 (5tn) + e g u (5»n) + • • • 

g 2 B g 20^* n ^ + e g 21^* n ^ + • * * * ( 4 *6) 

Again by applying the rules of differentiation, it can be shown that 

3Z _ 3T . 3K 
35 " 35 6 35 


a£z _ 3 2 T . 3 2 K 

U 2 ' 5F * sF 


3 2 Z _ 3 2 T 3 2 K 

3£3n ” 353n 6 353n 


(4.7) 


where 


Z = f. 


6i * 8<" 


T 




g 10 » &20 


and 


K = f ll* f 21* g u» g 21* respectively. 


Substituting the expressions of (4,6) and (4.7) into equations (4.5) and 
keeping terms only of 0(1) and 0(e), the resulting equations become 
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3 2 f 

IP 


10 


J 2 f 

♦ f m + + f,, + 2 


3 2 f, 


af. 


10 ' -5I 2 — Ti ii T <7&r r a n 


^ + o3 +2f 


if, 

20 Jf 


10 


3f 20 3g 10 3g 20 a2f ' 
+ 10 35 + 2g 20 3g + 2g io 35 * K 35^ 


+ 10*20 + g 10 g 20^ 8 0 


a 2 g 

K* 


10 


+ + eb 


a 2 g 


a2g 10 _ 3g l 


11 . _ . - *10 . - "•10 . "‘10 


3f, 


10 * eL W~- + g u + 2 mr + 0 35“ + 2g 


20 35 


3g 20 3f 

+ 2f 10 3 T" ' 2g 10 TT" ■ 2f 20 +K W 


20 “20 „ 3g 10 . „ 32g 10 


+ 2nw(f 10 g 20 - f 20 g 10 )] 8 0 


3 2 f a2f 

W~~ + 4f 20 + e[— g— y2. + 4f 21 + 2 


32f 20 , -**20 
3£3n a W~~ 


3g 3f g2f 

+ g 10 35 f 10 35^ + 4K 35^ + £ n(g 10 2 - f 10 2U = 0 


a 2 g<- 


3 2 g, 


HT- + 4g 20 + + 4g 21 + 2 7&rT + ° “ &2 ° 


32g 20 — 3g 2 

It 


3f 


3g, 


"8, 


10 _ f "*»10 A llV " *20 - 


3 2 g, 


10 — - f io — + 4K -3 jr- 8 0 


(4.8) 


By separating the terras of 0(1) and 0(e) in the equations of (4.8) and 
equating both sets of terms equal to zero, the resulting equations become 


32f]L0 + f 

gp— + f 


10 


* 0 
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35 2 " 

* 2f 20 

W~ 

***20 

w ~ 


* *10 " ° 


+4f 20 * 0 


+%0 9 0 


(.4. 9. a) 


* 2 'll + f s .^flO _ -iflO f fl0. 2f ^20 
f ll *353^ a 35 2f 2035 2f 103 5 


3«10 ®*20 3 2 fi 0 

” 2 * 20*35 2 *10 *3^“ “ K 7^ 2n <^ f lo f 20 + *10*20* 


52 ®u 

+ *11 * ‘ 

o 

3 2 *10 

- 3 *10 

35* 


353n 


+2 *10 

3 ^20 

3 *10 

- K ^ 

85 “ 

20 

35 

K 85^ 

3 2 f 2i 

t 4f n = 


3 2 f 

2o 

- » 2f 20 

3C Z 

-2 

353n 

0 35 

j2f 

" 4K W~ ’ 2 wn(g 10 2 ~ ; 

f io 2) 

3 2g 21 

x ii n * 


32g 20 

- 3 *20 

35 2 

+ 4g 2 i 


353n 

" ° 35 


2g *“- 2 f 1S2. 

^*20 35 "lO 35 


2n ^f 10 *20 ■ f 20*10* 


3 *10 - 3f io 

‘ *10 35 * r 10 35 


+ g 


3 f 

io IT 


10 


3g 


10 


3 2 g 


20 . 


+ f 10 35 " 4K 3 f 2 + n “> (f io*10* • 


35 


(4.9.b) 
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The equations of (4.9. a) are linear second-order differential equations. 
Therefore, it can be shown that assuming the appropriate form of a solu- 
tion, the results become 

f 10 8 A^(n) cos 5 + Bj^Cn) sin £ 

g 10 “ A 2* n * 008 5 + B 2* n * 8in 5 

f 20 * A 3^ 008 2 * + B 3^ n ^ sin 25 

*20 “ A 4^ n ^ cos + B 4^n) sin 2 5 . (4.10) 

Substituting (4.10) into (4.9.b) and using the multiple- angle formulas 
yields 


3 2 f 


11 


+ f., s -2 


35 


2 11 


dA x i _ i 

IT - 7 V + f A 1 A 3 + B lV 


-(B^ + A : A 3 ) + i(A 2 A 4 + B 2 B 4 ) - (A 2 A 4 + B^) 


+ nuC^A^g - A 3 B 1 ) + i(A 2 B 4 - B^))! sin 5 


’] 


-2 


dB l 1- 1 

ar + r B i ♦ f Vi - A iV + (a i b 3 - Vi> 


+ \ (A 4 B 2 - A 2 B 4 ) + (A 2 B 4 - B 2 A 4 ) - iKA 1 


+nu)[-i(A^A 3 + B^Bg) + y (A 2 A 4 + ® 2 B 4 )]J cos £ + . . . 


* 2 gll 

H z + g n " “ 2 


’ dn “ 2 ^ A 2 + 2^ A l A 4 + ®lV 


-( A ^ A ^ ♦ + ^ A 2 A 3 + B 2 ® 3 ^ “ 2 ^ A 2 A 3 * ® 2 ® 3 ^ 

- i KB 2 + nStjU^ - Ufa) - i(A 2 B 3 - A^jj sin l 

-{a^ * ¥ b 2 * ¥\h * Vi> ♦ <*1 S 4 - w 


-<*2 B 3 - B 2*3> - 1< A 3 B 2 - W ■ T “2 


nu rl 


C ^< A ^ A 4 ♦ “ ^ A 3^2 + ® 2 ® 3 ^ COS ^ + • • • 

m 


3 2 f 

21 

sV 2 + 4fzi * ’ 2 


+ i [i(B 2 2 - A 2 2 ) - ^2 - A^)] + 2K(-4B 3 ) 

K ca 2 b 2 - A iv] sin « -{4r * 


“3 1 - 

*drT " 2° (2A 3 ) 


+2K(-4A 3 ) ♦ ^»w[i<A 2 2 - B 2 2 ) - ^(A^-Bj 2 )] 
* 2 g 


cos 2£ + . . • 


H 


*21 

r + 4g 2i s - 2 


-*&T * r C_2A 4 3 ■ H< B 1 B 2 " A lV 


♦ ^ B 1 B 2 ' W 3 + 2KC ' 4 V ■ ^>t|<AiV A 2 B l )3 

f dB 

? t2 V - ^ WVl ’^ WW 3 


sin 2£ 


43 


♦ 2Kt-*A„J - inSiUjAj 


»1« 2 )J 


CO* 2ft*.. 


C4.ll) 


where ♦ . . . indicat** terms multiplied by sine* and cosines of integral 
multiples of g other than those shown. The particular solutions corre- 
sponding to the terms shown on the "ight-hand sides of (4.11) will contain 
terms proportional to g sin ng or g cos ng [n * 1 for (4. 11. a, b), n ■ 2 
for ( 4 . 11 . c , d)]. Thus, the second approximation would be unbounded for 
large g while the first approximation is bounded for all £. These 
unbounded terras are called singular terms. The terms on the right-hand 
sides of (u.ii) indicated by + ... do not lead to singular terms. 

The idea of a perturbation solution is that higher order terms in 
the series solution represent small corrections to this first term to 
obtain a uniformly valid expansion. The presence of this singularity 
causes this fundamental idea to be violated. Therefore, since the expres- 
sions of n dependency are independent of the variable causing the singu- 
larity, the n-dependent expressions can be set individually equal to zero 
to avoid this problem. Therefore, from equations (4.11), the resulting 
equations, which are eight ordinary first-order differential equations 
having n dependency, become 


dA 


i— i i 

IT + r A l + T*l * ^ A l A 3 + B 1 B 3 * A 2 A 4 + B 2 V 


♦ ■jnwCBjAg - AjBj + BjA 4 - A 2 B 4 3 = 0 


dB l 1- 1 1 

dT" * 2° B 1 - 2**! + ? A l B 3 - B 1 A 3 “ A 4 B 2 + W 


t t ♦ BjB 4 ] * 0 

JT ♦ ♦ jKBj ♦ ^V„ ♦ B 1 B 4 - AjAj - B 2 B,3 

+ -|n«t-A^B 4 + A^B^ t AjBj - AjBj] * 0 

^®2 1—1 1 

♦ jO» 2 - -jKA 2 ♦ ^*A 4 B x + B 4 A x - A 2 A 3 t A^ 

♦ ■jftutAj^ ♦ B^S 4 - AjAj - BjBj} « 0 

dA« 1 \ 

ST * ?“*» * *»> ♦ T t» 2 2 - >I 2 ♦ v - V3 

t ^unCAjBj. - A 2 B 2 3 e o 
dB, , , 

JH- * 1^3 - 4KA 3 * ? A 2 B 2 • W 

♦ -JgSlCAji - Bj2 - Aj! ♦ B^] ■ 0 

1— 1 

ZT * 2° A 4 + 4KB 4 + 4 {B 1 B 2 " A lV 

+ ^wnCA^Bj + AjB^] = 0 
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dB„ 


3»4 - 


4KA., 


frVi 


A 2 B 1 3 


- ^UnCA.Aj - B^ap * 0 


(4.12) 


Since aquations (4.12) arc first-order nonlinear ordinary differential 
equations, tha fourth-order Runge-Kutta program, previously developed, 
can be used to solve for the modal amplitude coefficients. By finding 
these coefficients for various points in time, a relation between the 
results of equation (3.21) and equation (4.12) can be observed to the 
approximation of order e. 

Solving equations (4.12) for the highest derivative (first order 
ir .his case) and substituting n s et, the governing equations for the 
Funge-Kutta program become 


1— 1 1 

dt~ “ C 2° A 1 " 2^1 " 2* tA l A 3 * B 1 B 3 + A 2 A 4 + B 2®4^ 


1 


- fn^( Bl A 3 - Al 3 3 ♦ B 2 A 4 - A^)] 


dB. 


dt ' c 2 c3 1 + 7^1 " 2^ A 1 B 3 “ B 1 A 3 “ A 4 B 2 + A 2 B 4* 


- ♦ B 1 B 3 + AjA 4 + B 2 B 4 )] 


^ A "5 1_, 1 1 

dt” ' c 2° h 2 ~ r 3 2 ‘ 2* A 1 A 4 + B 1 B 4 - A 2 A 3 “ B 2 B 3* 


- ■jnu)(A 4 B^ - A^B 4 + AjBg - AgBjJJ 
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- ■|wn( J A 1 A 4 t BjB 4 - AjA 2 - BjBgJJ 


dt = eC “ 2° A 3 " 4KB 3 “ 8 <jA 2 2 “ B 2 2 + B l 2 ” A l 2) 


- inuCA^ - A 2 B 2 )] 


dS 3 x _ ! 

d T = eC " 7*3 + 4KA 3 “ ? A 2 B 2 - A lV 


j|iwo(A 2 2 - B 2 2 - A^ 2 + B^ 2 )] 


<*\ x _ ! 

dt = e ^“ 2° A 4 “ 4 ‘ ?:1 4 " 4^ B 1 B 2 ” A 1 A 2^ 


- •|nu)(A 1 B 2 + AgB^)] 


!_ i 

IT = e[ - H + 4KA 4 + ?Vl + A 2 B 1 ) 


+ ^unCA^ - B^)] . 


(4.13) 


It is often convenient to express the equations for A^ and B^ in 
terms of amplitudes, C^, and phase angles, which are also functions 

of the slow time variable n. Mathematically, we can express the relation- 
ships between the quantities as 


A. = C. cos 4. 
l i 


(4. 14. a) 


B. = C. sin <4. 
l i r i 


(4.14.b) 
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dA, dC, 


d*. 


dT a dT cos *i " c i dT sin *i 


C4.14.C) 


dB. dC. 

3T s dT sin *i + c i dT cos *i 


C4.14.d) 


where i = 1, 2, 3, and 4 for each of the equations above. Substituting 
the expressions of (4.14) into the first two equations of (4.12), the 
resulting equations become 


cos - C L Jp sin ^ oC x cos ^ + ^sin ^ 


+ itC^Cg cos 4> 1 cos sin ^ s m + 3 


+C 2 C 4 cos cos * 4 + C 2 C 4 sin * 2 sin * 4 ] + ^ c i c 3 
sin cos * 3 - C^g cos ^ sin t 3 + C 2 C 4 cos $ 4 sin <fr 2 

- C 2 C 4 cos d> 2 sin <j) 4 ] = 0 

dC d$. l— 1 

— - sin + C x cos sin <j> 1 - j^Cj. cos 

•+ cos <^ 1 sin 4 > 3 - cos <f> 3 sin ^ 

- C 2 C 4 cos $ 4 sin $ 2 + cos $ 2 sin * 4 3 + 

cos (j >1 cos <t> 3 + C^g sin ^ sin $ 3 + C 2 C 4 cos $ 2 cos $ 4 



i 


+ C 2 C 4 sin <fr 2 sin ^ = 0 
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(4.15) 


Multiplying the first equation by cos and the second equation by sin^^, 
adding the two expressions together, and vising appropriate multiple-angle 
identities from trigonometry, the resulting equation for becomes 

dC l 1 - 1 

dn“ + 2° C 1 + 2 fC l C 3 Ccos(24, l “ *3 )] 

+ C 2 C 4 [cos ($ 2 - + 4 + + ^nwtC^Cj 

sin(2<j > 1 - 4> 3 ) + C 2 C 4 sin (<|>2 - 4 > 4 + = 0 . (4.16) 

Similarly, multiplying the first equation of (4.15) by -sin and the 
second equation by cos <^, adding the two expressions together, and using 
appropriate multiple-angle identities for trigonometry, the resulting 
equation for becomes 


iiT ' ¥ ' ¥ C 3 sin(2 *l ' V + ~ si n (* 2 "V*l )] 

c c 

+ -|nu)[C 3 cos(24> 1 -<(> 3 ) + 4 cos(4> 1 +<j>2-4' 4 )] = 0 . 


(4.17) 


Using these procedures discussed above, equations for Cj, ^ 3 * $ 3 * C 4 » 
and can be derived. Thus, these transformed equations are 


dC 2 1 - 1 

dT + 2° C 2 + cos <*r* 4 + *2 ) 




+ •intjjCC 1 C 4 s in( 4> 1 -4> 4 +4> 2 ) - C 2 C 3 sin( 24 > 2 - < l» 3 )^ = 0 
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^$2 1 1 ^1^4 

dT ' 2 K - ■i»C4 1 -# 4 +# 2 ) - C 3 .in(2# a -# s U 

c c 

+ |nu[-~- cos(^ 1 -# 4 +^ 2 ) + C 3 cos(2$ 2 -$ 3 )J = 0 
dC 3 1- 1 

d n" + 2° C 3 + 8^C 2 2 c0s( 2 *2"*3 ) ‘ C 1 2 cos(2^ 1 -^g)J 
- •j|n(u[C 2 2 sin(2<(> 2 -$ 3 ) - C 1 2 sin(2(j> 1 -$ 3 )] = 0 
d$ C 2 C 2 

dT ‘ 4K + ®inC2# 2 -* 3 ) ~ cT" sin(2* 1 -* 3 )] 

0 3 

c 2 c 2 

+ cos C 2< ^2“*3^ “ r~ cos^S^-^g)] = o 

3 3 

dC 

dn~ + ^° C 4 " 'JtC 1 C 2 cos(^ 1 +^ 2 -^ 4 )] + inSCCjCj 
sin^ + <f> 2 - i^ 4 )3 = 0 


<*$ 4 

dVr 


1 ^1^2 
- 4K + ^"5~ 


c c 

sin(^+^ 2 -^ 4 >] - 


cos^ + <J> 2 - <j> 4 >] = 0 . 


(4.18) 


Equations (4.16), (4.17), and (4.18) are the general combustion 
equations in terms of amplitudes and phase angles. From this point, 
special cases can be investigated isolating certain conditions and closed- 
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form solutions can be obtained for these cases. It is convenient to do 
this in order to check the closed-form results of the special cases with 
the results from the general equations (4.16), (4.17), and (4.18) when 
the same conditions are imposed. 

The first case to be evaluated is the case for standing waves 
with no combustion effects. To simulate standing wave effect, set the 
amplitudes C 2 and and phase angles <J> 2 and equal to zero. This 
automatically satisfies four of the eight equations (4.18). To achieve 
the no-combustion effect, set the interaction index, n, equal to zero. 
Also, set the correction variable, K, equal to zero since the effect of 
K will be investigated separately at a later time. Imposing these con- 
ditions, the governing equations reduce to 

dC l 1 - 1 

dn + J 00 1 + 2 C 1 C 3 = 0 (4. 19. a) 


~ ^sinU^-^n = 0 (4.19. b) 


3 1-1 

d “ + 2 ° C 3 " S’ c i 2cos ( 2 $i"$3) = 0 
dl ^3 1 ^ 1 ^ 

= °* 


(4.19.c) 


(n.l9.d) 


The initial conditions imposed for this case are 


^( 0 ) = 1 
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c 3 (o) * 0 

♦3/0) = 4 > 10 

*3 (0) = *30 . 

To attempt a closed-form solution , let 



Substituting these expressions into equations (4. 19. a) and (4 
dividing through by e~^° n , the resulting equations become 


fi ♦ i =os(2 h - - ° 

Multiplying equation (4.2 2.a) by 1/4 and equation (4.2 2.b) by 
adding the two equations, terras containing the cos(2<ji^ - <j>^ )e 
eliminated. In doing so, the result becomes 


(4.20) 

(4. 21. a) 
(4.21.b) 

(4.21. c) 

(4.21.d) 
19. c) and 

(4. 22. a) 

(4.22.b) 

F^/? 1 and 
■%5n 


are 
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cn 


F 3 dF 3 


T 1 da 


0 . 


(4.23) 


Multiplying through equation (4.23) by F^ gives 

J S' CF i 2 * “V 3 * 0 • < 4 - 21 *) 

Integrating with respect to n then dividing by 1/2, the resulting equa- 
tion becomes 


F l 2 + 4F 3 2 = D 1 (4.25) 

where is a constant of integration. This constant depends upon the 
initial condition imposed on the problem. From the initial conditions 
given in (4.20) and using the transformation (4. 21. a) and (4.21.b), it can 
be shown that F^(0) = 1 and Fg(0) = 0. Therefore, equals to 1. Thus, 
equation (4.25) becomes 

Fj* = 1 - 4F 3 2 > (4.26) 

Taking equation (4.26) and substituting into equation (4.22.b), then 
separating variables, the resulting equation becomes 

3_ = 1 “%on cos ( 2 <j) )dri (4.27) 

[1-4F 3 2 ] 8 13* 

Letting 2^ - $ 3 = -tn, which satisfies equations (4.19.b, d), yields 

l 

cos( 24 >^ - <j> 3 ) = (-1) where £ = 0, 1,2,3. . . Substituting this expres- 
sion and integrating the above equation, the resulting equation becomes 
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y tanh -1 2F 3 » i [- | e'^+DjJC-l)* (4.28) 

where is a constant of integration. Using the initial condition Fg(0)« 
0, then, it can be shown that = 2/3. Substituting and taking the 
hyperbolic tangent of both sides of equation (4.28), the result becomes 

F 3 =Vtanhti-(-l)^(l-e"^ n )] , (4.29) 

Substituting this expression into equation (4.26) and simplifying, the 
resulting equation becomes 

F = sechC^- (l-e’ %an )3 . (4.30) 

1 2o 

Substituting equations (4.29) and (4.30) into equations (4. 21. a) and 
(4. 21.b) , and substituting n = et and u = 5e, the resulting closed-form 
solution for wave amplitudes and are 

C = e"^ t (sech[ E ^ )£ ( l-e~^ Mt ) ] ) (4. 31. a) 

2w 

C = . (4.31.b) 

To find expressions for and <|> 3 , substitute the relation that 2<j>^-4> 3 =£iT 
into equations (4.19.b) and (4.19.d) and integrate and evaluate the con- 
stants of integration with the initial conditions; the results are 


♦l = *10 


*s * * tha -* 1 
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C4.32) 

u 

where is a constant and is tv radians out of phase with 2$^. It 
can be seen that a special set of initial conditions is necessary to be 
consistent with this solution. A representative set is ^ = $g Q *• 0 
which corresponds to t * 0. 

Inspection of equations (4.31) reveals that the magnitude of 
continually decreases with time while the magnitude of Cg first increases 
and then decreases. An interesting special case of equations (4.31) 
occurs in the absence of steady-state combustion (it a 0). The results of 
this case are 


C L = sech[ ( ~ 1 ^ - ”] 

t 

C 3 = y tanh C— E — ] ^ (4.33) 

These results show that a disturbance in the form jf the first mode is 
transferred to the second mode as time increases. It is thought that this 
indicates the beginning of the steepening that leads to the formation of a 
shock wave. It can be seen that the presence of damping, in the form of 
steady-state combustion, inhibits this process. 

The second case to be investigated is that of standing waves with 
gas-dynamic nonlinearities neglected. To simulate the standing wave 
effect, let the amplitudes Cg and and the phase angles $g an< * $4 ec l ua ^ 
zero. Again, this automatically satisfies four of the eight equations of 
(4.18). To achieve omission of gas-dynamic nonlinearities, let i = 0. 
Also, let the correction variable, K, be equal to zero for simplicity. 


In doing so, the resulting equations, based upon equation (.4.18), 
become 
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oC^ + y nwtC^Cj sin^ - * 0 

nuCC 3 cos(2^ - $ 3 )D = 0 

oC 3 - -j|- nw[-C^ 2 sin(2^ - $ g )] = 0 
d 4> 3 i _ C 2 

^ n«C-^ cos(2^ - $ 3 )] * 0 . (4.34) 

3 

The initial conditions imposed for this case are 
C^O) = 1 

C 3 (0) = 0 
♦l< 0) * *10 

♦,(0> ■ *30 . <4 - n> 

Let 2^ - $ g * (2^ + l)ir/2, t * 0, 1, 2 ... . This implies that sin(2^ 1 
! 

- $ g ) * (-1) and 008(2^^ - #3) = 0. Substituting into (4.34) and solving 
in the manner indicated previously one obtains expressions for the ampli- 
tudes for and C g which are 

s e "^“ t (secC-^ ne(-l)^(l-e’^ t )]} 


flC l . 1 
TT + 2 


d *l 1 

W+ 2 


<ic. 


(4. 36. a) 
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C s ncC-l^a-e^JJ} (4.36.b) 

u 3 T7T H 

♦ x * ^ 10 (4.36. c) 

♦3 ■ s *io - (£ r i >' . (4 - 36 - d > 

where is constant and $ 3 is (2fc+l)v/2 radians out of phase with 2^. 

As in the previous solution, special initial conditions are required to 
produce this solution. A representative set is ^ * 0, $ 3Q * -v/2, which 
corresponds to t * 0. 

The secant and tangent both become infinite when their arguments 

take on the value *v/2. In (4. 36. a, b), the arguments of these functions 

3/2 

start at zero at t * 0 and have a maximum absolute value at ne/2 . 

3/2 

Thus, if ne/2 < v/2, the tangent and secant never become infinite and 

and C 3 eventually decay to zero due to the influence of the exponential 

3/2 

function. This is a stable situation. If, on the other hand, ne/2 > 
v/2, the tangent and secant become infinite at t w 8 (2/w)|£n[l-2^v/(ne)]| 
causing and C 3 to become infinite. This is an unstable situation. 

Thus, the boundary between stable and unstable behavior is indicated by 
the equation 


ne/2 3/2 8 v/2 . (4.37) 

The stability equation in the n-e plane has the form 

n 8 2^v/e 8 4.442/e. (4.38) 


This has the form of a rectangular hyperbola and is independent of w. 
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For the cam of traveling vay«i, it it sort convenient to work 
with tha general perturbation equations expressed in nodal amplitudes 
in terms of the real time variables, equation (4.13). To simulate the 
effect of spinning or traveling waves, let the following modal ampli- 
tudes be equal. These relations are 


B 


1 



»3 ■ A . < “- 39) 

It can be shown that substituting the relations (4.39) into equation 
(4.10), expressing the results in terms of the real time variables, sub- 
stituting these expressions into equation (3.18), and using appropriate 
multiple-angle formulas leads to 

♦(0,t) * Ajcos(t-e) - A 2 sin(t-©) ♦ AgCos 2(t-9) 


•A^sin 2(t-0) + . . • . 


(4.40) 


which has the form of a sum of traveling waves. Substituting the expres- 
sions in (4.39) into equations (4.13), these eight equations reduce to 
four pairs of identical equations. The four independent e nations listed 
below are 


dAi 

IT 


c[- y 0 A 1 +iKA 2 -i(A 1 A 3 +A 2 A 4 )-nJ(A 1 A 4 -A 2 A 3 )] 
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3“ * c[- 7 ® A 2 “ ? ^l“ i ^ A l A 4" A 2 A 3^ +n “^ A 2 A 4 +A i A 3^ 

dA, , . , , 

■ c[- -j aA 3 +4KA 4 + -jjKA^ 2 ^ 2 )-*— nuCA^A^JJ 

flT * e[- j 5A 4 -4KA3+ j KA^)- ^*to(A 1 2 -A 2 2 )] . (4.41) 

By making the substitution, we have reduced to a system of four equations 
and four unknowns. By solving for the modal amplitudes A^, the modal 
amplitudes B^ are readily computed by using the relations of (4.39) to 
determine the entire nature of the wave form. 

Tor the case of traveling waves omitting gas-dynamic nonlinearities, 
let the amplitudes A^ and A 3 equal zero. T- set 1, the gas-dynamic 
index, equal to zero. Again, for simplicity, let the correction variable, 
K, controlling physical chamber configurations, be zero. In doing so, in 
terms of the transformation variable, n, the resulting equations become 



- nu[A 2 A 4 ] * 0 




0 


(4.42) 


which Is a system of two equations and two unknown modal amplitudes. To 
find an exact closed-form solution to these equations, let 


A 2 



A 


4 
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a, fl ♦.-*>£ 

dA. dF. 

__ . | 5)F 2 ♦ ^ . (4.43) 

Using these transformations, the procedure for solution is exactly the 
same as for the standing wave case for both no combustion and no gas 
dynamics. The initial conditions for this case are 

A 2 (0) * 1 

A 4 (0) * 0 , (4.44) 

Substituting the expressions of (4.42) into (4.41), the resulting equa- 
tions are 


ar - .-v 5 " . o 


with initial conditions 


r i<°> • 1 

Fj(0) * 0 


Solving these equations in the manner outlined in the standing wave solu- 
tions, the results are 
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f, = s«c44 
a 


T. = . 


2 2/2 2 a 


(4.46) 


Expressing the results of (4.45) in terms of modal amplitudes by substi- 
tuting into (4.42), the resulting equ^ions become 


A 2 = si. 


A 4 


2/2 * a 


(4.47) 


The results for traveling waves (4.47) are quite similar to the results 
for standing waves (4.36) for the case of no gas-dynamic nonlinearities. 
The same behavior can be expected as was discussed in the standing wave 
case about the nature of oscillation of the modal amplitudes. The only 
significant difference is the value to deterr-e the boundary of stability 
for the interaction index governing the combustion terms. The stability 
condition for traveling waves is 

/? r _ ir 

” - n£ “ t' 

2 2 • (4.48) 

Thus, the equation of the stability boundary in the n-e plane is 

_ it 2.22 

11 " " e • t 4 - 4 9> 


Comparing equation (4.49) to (4.38) shows that the stability boundary for 
the interaction index is half as great for the traveling wave case as for 
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the standing wave case for any e. This will be yerified in a later pre- 
sentation of results of yarious numerical cases. 

For the case of traveling waves with no combustion, let the ampli- 
tudes Aj and equal to zero. Then set n, the interaction index, equal 
to zero, and, again, let the correction .variable K equal to zero. Sub- 
stituting into equations (4.40) and transforming into variable n» the . 
results are 


^1 1 - 

dT + 2 ffA l + A 1 A 3 = ° 


f^3 . 1 - . 1 * 2-0 

dn * 2 0 A 3 " 4 A 1 ‘ 0 


with initial conditions 


A x (0) = 1 

a 3 (o) = 0 


(4.50) 


which again is a system of two equations and two unknown modal amplitudes 
To find an exact closed-form solution to these equations, use similar 
transformations as shown in (4.42). In doing so, and simplifying, the 




results are 
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with initial conditions 


F^O) = 1 

f 2 (o) = 0 . 

Solving these equations in the same manner as before, the results are 

F^ = sech[l/o(l-e""^ n )] 

F 2 = | tanh[l/5(l-e' %5n )3 . (4.52) 

Again, expressing the results of (4.51) in terms of the modal amplitudes 
of the form of equation (4.43), the resulting equations become 

A^ = e”^° n sech[l/a(l-e ^ an )3 

A = tanh[l/o(l-e~^ on )] . (4.53) 

The results for the traveling waves (4.52) are similar to the results 
for standing waves (4.31) for the case of no combustion. A disturbance 
initially having the form of the first mode eventually is transformed into 
one having the form of the second mode. To compare these results for 
standing waves and traveling waves to the general perturbation equations, 
two computer programs were written (Appendices D and E) which numerically 
evaluate the modal amplitudes of various conditions for standing and 


traveling waves. 
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One last special case is an investigation of the effect of the 
correction variable K. In the special cases previously discussed, the 
correction variable K was set equal to zero. But, in this discussion, 
the correction variable K will be of primary importance in the equations. 
To start this analysis, refer to equations (3.21). Based upon these 
equations, impose the following conditions. First, neglect combustion 
effects (i.e., n = 0). Then, let us consider only the case of standing 
waves (i.e., - 0). Finally, let us neglect the steady state 

burning rate (i.e., o = 0) and assume that the terms multiplied by eK 
are larger than those multiplied by c above. This can be accomplished 
by writing 

K. = eK (4.54) 

and treating as a quantity of 0(1). Imposing the above conditions 
and substituting equation (4.54) into the equations (3.21), the result- 
ing equations become 


d 2f l 

df l 

df 2 


[1+K.] - + f, + 2e[f, 

1 dt* 1 2 

dT + 

f l dT"- 1 = 0 

(4. 55. a) 

d 2 f „ 

df i 



[1+4K.] - + 4f- - ef, 

1 dt 2 21 

X _ 

dt 

0 

(4.55.b) 


with initial conditions 


f x (0) = 1 
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df. 


ar< 0) - 0 

* 2 ( 0 ) = 0 


df, 


dt 


:o> = 0 . 


First, assume a straightforward perturbation solution similar to the 
equations (4.6) except the functions are dependent upon the real time t. 
Substituting these assumed solutions into the equations and initial con- 
ditions of (4.55) and keeping termc of 0(1) and 0(e), the separated 
equations become 


d 2 f 


10 


dt2 + TI^T f io 


= o 


(4. 56. a) 


d2f 20 4 

dt 2 n?4K]T *20 ‘ 


(4.56.b) 


d 2 f 


11 + (TT^)f, 


df 




10 


df, 


dt 2 1+K^ 11 (H-K ) u 20 dt 


- f 


20 - 


10 dt 


(4.56.c) 


d 2 f 


df 


21 . /_4 . , . 1 r r ”10- . 

d 2 '1+41^' 21 " TI+4K^T- 10 dt J 


(4.56.d) 


with initial conditions 


f io«» - 1 


f u (0) = 0 
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0 


-~=< 0 ) = 0 


f 20 (0) =0 



Co) * 0 


dt 


= 0 


df 21 (0) 

dt 


= 0 . 


The first-order equations (4. 56. a and b) can be solved by assuming the 
usual assumed solution for linear differential equations. Doing this 
and applying the appropriate initial conditions, the results for the 
first-order terms are 




(4.57) 


Substituting (4.57) into the right-hand side of (4.56.c) the equation 
becomes a homogeneous linear differential equation. Solving in the 
usual manner and applying the appropriate initial conditions 

f u = 0 . (4.58) 


Substituting (4.57) into the right-hand side of equation (4.56.d), the 
resulting equation becomes a linear differential equation with a particu- 
lar solution. By assuming an appropriate homogeneous and particular 
solution and evaluating the constants using the appropriate initial condi- 
tions , the result becomes 


2 


2 


(4.59) 



-i'i+4K 1 

24K X 


sin- 


t + 



sin 


Therefore, substituting equations (4.57), C4.58), and (4.59) into the 
assumed perturbation solution and letting = eK, the resulting equations 
become 


f^ = cos 


(4. 60. a) 


f 2 = 


24K 


tr 


l+4Ke 

1+Ke 


sin 


/l+4ek 


t - sin 


/I+kF 


t] . 


(4.60.b) 


Recall that in the two-variable perturbation method, and expressed 
in terms of the perturbation variables were 

f^ = A^(n) cos £ + B 1 (n) sin £ (4. 61. a) 

f 2 = Ag(n) cos 2£ + B 3 (n) sin 2£ . (4.61.b) 

By transforming equation (4. 60. a) into perturbation variables and expand- 
ing the argument of the cosine function by the Taylor series and using 
appropriate sum and difference trigonometric identities f^ can be 
expressed as 

f^ = cos ~ Kn cos £ + sin y Kn sin £ . (4.62) 

Therefore, comparing this to equation (4. 61. a), the functions and B 1 


must be 
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A^(n) * cos ^ Kn 

B^n) = sin j Kn . (4.63) 

By similar procedure, it can be shown that evaluating equation (4.60.b) 
and comparing it to equation (4.61.b), the results are 

Aj( n) = Tj^Csin 4Kn - sin Kn3 

B 3 (n) = 2 ^ cos Kn " cos 4Kr ^ • (4.64) 

To show the validity of equations (4.63) and (4.64), the problem is now 
solved using equations (4.12) which are derived from equations (3.21) by 
the use of the two-variable perturbation method. To reproduce the condi- 
tions imposed on the problem just discussed, let there be no combustion 
(i.e., n = 0), let there be no steady-state burning rate (i.e., c = 0), 
and let there be only standing waves existing (i.e., k^ = k^ = B 2 = B^ = 
0). Imposing these conditions on equations (4.12), the resulting equa- 




dA, 

~ + 4KB. + 
dn 3 


+ 2 ^ 1^3 

* K B 3 

I^ 2 - 


♦ B 1 E 3 ] * 0 


B lV ' 



0 


tions become 
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“a 1 

ZT- 4KA 3-7 A 1 B 1*°- 


(4.65) 


In the previous solution it was assumed that the frequency correction 
terms were larger than the gas-dynamic .nonlinearities. To be consistent 
with this assumption the following procedure is used. By a change of 
variable n = ;/K, equation (4.65) can be rewritten as 


a*i i i 

dc + 2 B 1 + 2K CA 1 A 3 + B 1 B 3 ] 


= 0 


dB 

7 A 1 + 2K CA 1 B 3 " B 1 A 3 ] 


= 0 


dA, . 

ST + 4B 3 + 8K CB 1 2 ~ A l 2] = 0 


^3 ! 

d; " 4A 3 ' 4K A l B l = 0 • 


(4.66) 


Assuming a straightforward expansion of the form 


A 1 = A 10 + K A 11 + * * * 


B 1 5 B 10 + K B 11 + * * ' 


A 3 = A 30 + K A 31 + * * * 


B 3 = B 30 + K B 31 + * ‘ * 


(4.67) 
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then substituting these expressions into the equations 04.66) and keeping 
terms of OCl) and 0(1/K), the resulting separated equations become 




(4. 68. a) 


dB 

10 1 . . 
d; " 2 A 10 " 0 


(4.68,b) 


^30 

"dc + 4B 30 3 0 


(4.68.C) 


^ B 30 

TT- 4A 3o s 0 


(4.68. d) 


1 

d; + 2 B 11 


2^ A 10 A 30 + B 10 B 30 ] 


(4.68.e) 


dB ll 1 
d; “ 2 A 11 


" 2* A 10 B 30 " B 10 A 30 ] 


(4.68.f ) 


dA 


dT + 4B 31 ■ - ^ B 10 2 - V ] 


(4.68.g) 


dB 


d? " 4A 31 = 4 CA 10 B 10- 1 


(4.68.h) 


with the initial conditions 


V°> * 1 B 10 (0) = 0 

A n (0) = 0 B 11 (0) = 0 
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W 0) * 0 B 3Q C0) 8 0 

A 31 (0) « 0 B 31 C0) « 0 . 


Since the first-order equations are coupled, differentiate equations 
(4, 68. a and c) once with respect to z then substitute equations (4.68.b 
and d) into these equations resulting in 


d * A 10 i 

dz 2 + 7 Al ° ‘ ° 


d 2 A 


30 


dc 2 


♦ 16A 30 8 0 


(4.69) 


As can be seen, equations (4.69) are linear differential equations 
which can be evaluated by the usual manner. In doing so and applying 
the appropriate initial conditions, the resulting first-order modal ampli- 
tudes are 


A lo = ’ CM 

A 30 • 0 . (4.70) 

Knowing values for A^g and Agg, substitute these values into equations 
(4.68.b and d) and apply appropriate initial conditions. The results 
become 

B 1q = sin jz - sin -jKn 


(4.71) 
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Substituting ths results of (4.70) Mid (4.71) into ths right-hand side of 
aquations (4.68.e-h), ths resulting equations become 


**11 . 1 B 
IT" T B 11 


(4. 72. a) 


TT 7 ii 


(4.72.b) 


^31 1 

TT + 4B 31 8 8 008 5 


(4.72.C) 



4A 31 * ? Sin 


(4.72.d) 


Since equations (4.72.C and d) are coupled, differentiate both equations 
once with respect to c and substituting equations (4.72.C and d) into the 
appropriate terms of the new set of equations, the resulting equations are 


d 2 A 


31 


dC 2 


♦ 16A 31 - 


| sin C 


d 2 B 


31 


de 2 


+ 16B 


31 


5 

? cos 


(4.73) 


Equations (4.73) are a set of linear differential equations with homo- 
geneous and particular solutions. Solving these equations in the usual 
manner and using the appropriate Initial conditions, the resulting modal 
amplitudes are 


A 31 8 2^ " sine) 8 -^<sin4Kn - sinKn) 


24' 
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®31 * 25 <c08C “ cos4 ^ ■ ^coiKn - cos4Kn' , (4.74) 

Zn a similar manner, the results for the modal amplitudes and can 
be determined to be 

A n . 0 

B u * 0 (4.75) 

evaluated with the appropriate initial conditions. Therefore, substitut- 
ing the results of (4.70), (4.71), (4.73), and (4.74) into the assumed 
perturbation solution of (4.67), the resulting modal amplitudes became 

* cos ^ Kn t # • • 

= sin Kn + . . . 

A 3 * 2^ sin 4Kn ~ sin + * • • 

Ej - -j^sin 4Kn - sin Kn) + . . . . (4.76) 

It can be seen that equations (4.76) are identical to equations (4.63) and 
(4.64). This indicates that the two-variable method produces the correct 
solution. Equations (4.60) indicate that the presence of K changes the 
frequency of each of the first two acoustic modes and further renders the 
ratio of the second frequency to the first a non-integer number in general. 


m 
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* 


* 


Equations (4.76) show how this offset manifests itself in the two-variable 
perturbation solution. 

These results can be used in another way. If the nonlinear terms 
are neglected in (4.S5.a) ( the results are 


(1+K X ) 


d*f 

i+ * 

dt 2 1 


0 


(1+4^) 



♦ 4f , 



0 


f x (0) * 1, 


df 1 (0) 

dt 


0, 


df-(0) 

f 2 (0) = 0, — * 0 . (4.77) 


It can be easily shown that equations (4.60) constitute the exact solution 
of equation (4.77). If the corresponding terms are neglected in equations 
(4.65), the results are 


dAi ^ 

dT + Pi ! 0 


dB i i 

dT-i^i 8 0 


dn 


+ 4KB g + •|<B 1 2 - A, 2 ) * 0 


d ®3 1 

dT - 4KA 3 " 4 A 1 B 1 = 0 


(4.78) 


where 
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A x CO) = 1 
B^(0) = 0 

a 3 (o) .= 0 
b 3 (o) = 0 . 

It can be shown that equations (4.76) are the exact solution of equation 
(4.78), These facts were used to check the accuracy of the computer pro- 
grams to be discussed later. 

In the remainder of this thesis, a comparison of the magnitudes 
of the modal amplitudes will be represented in graphical and tabular 
form. Under a given set of conditions, the acoustic modal amplitude pro- 
gram, the general perturbation program, and the analytical cases that 
were programmed will be used and results compared. Varying certain con- 
ditions will show their effect on the changes in magnitude of the modal 
amplitudes through a set range of time which is related to maintaining 
stability. From these various cases, it will be determined which param- 
eters and conditions have the greatest effect in changing modal ampli- 
tudes and which in turn affect the stability criteria for combustion 
by the methods discussed above. 


Chapter 5 


DISCUSSION AND PRESENTATION OF RESULTS 

In this chapter, results are presented both in graphical and 
tabular form which are representative of the results generated by the 
programs listed in the Appendices B through E. From these representative 
sets of results, basic observations will be made to observe which 
parameters or conditions have the greatest effects on the problems of 
stability. 

In Figures 3 and 4, modal amplitudes F^ and F 2 are graphically 
represented versus time for a stable standing wave case. For these 
figures, F^O) = 0, F^CO) = 1, F 2 <0) = 0, F 2 '(0) = 0, G 1 (0) = 0, G 1 '(0) = 
0, G 2 (0) = 0, G 2 '(0) = 0, n = 35, i = 1, K = 0, e = 0.1 and u = 0.1. 

The step size used was 0.1. Experimentation showed that this was a small 
enough step size to produce accurate results and was used throughout. 

From these figures, one notices that both the first and second order modal 
amplitudes decrease in amplitude with increasing time. Also, F 2 , the 
second order modal amplitude, tends to oscillate at twice the frequency of 
F^. These figures are based upon one set of parametric values; however, 
these figures represent qualitatively the results obtained using a wide 
variety of initial conditions and parametric values. In Figures 5 through 
8, modal amplitudes F 1 , F 2 , G^, and G 2 are grap’ ally represented versus 
time for a stable traveling wave case. For these figures, F^(0) = 0, 
F 1 '(0) = -1, F 2 (0) = 0, G 1 (0) = 1, G 1 '(0) = 0, G 2 (0) = 0, G 2 '(0) = 0, 
n = 15, i = 1, K= 0, u> = 0.1 and e = 0.1. The general shape of the 
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Modal Amplitude F„ vs time t for Traveling Waves - Stable Case 
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curves and the relative frequencies of oscillation are qualitatively 
similar to che stable standing wave case. 

In Figures 9 and 10, modal amplitudes F ^ and F 2 are graphically 
represented versus time for an unstable standing wave case with the 
same conditions as the stable case except that n = 50. As can be seen, 
the maximum amplitude of F^ starts to decrease then increase dramatically 
for increasing time. The maximum amplitude of F 2 increases continuously. 

In Figures 11 through 14, modal amplitudes Fj, F 2 , and G 2 are 
represented versus time for an unstable traveling wave case. Again, the 
conditions are the same as for the stable traveling wave case except that 
n = 30. Drastic increases in amplitudes are observed for all the modal 
amplitudes shown as time increases. The behavior is similar to the 
unstable standing wave case. The period of time for traveling waves to 
become unstable is about one-half the period of time for standing waves 
to become unstable. Thus, it seems that traveling waves are less 
stable than are standing waves. 

In Tables 1 and 2, a comparison of results is presented for modal 
amplitudes F^ and Fj for a stable standing wave case. For these cases, 
F 1 (0) = 0, F^iO) = 1, F 2 (0) = 0, F 2 '(0) = 0, G 1 (0) = 0, G 1 '(0) = 0, 

G 2 (Q) = 0, G 2 (0) = 0, n = 60, e = 0.1, and oj = 0.1. These tables 
quantitatively show the effect of neglecting gas dynamic non-linearities 
on the accuracy of the computations. Also, a comparison can be made 
between the exact solution method (Appendix B program) and the perturbation 
solution method (Appendix C program). From Table 1, one can observe that 
the effect of neglecting gas-dynamic nonlinearities is small where 
quantitatively comparing values of the modal amplitude F^. Even though, 
quantitatively, the values for the exact solutions and perturbation 




Amplitude F„ vs time t for Standing Waves - Unstable Case 







Amplitude F vs time t for Traveling Waves - Unstable Case 







Figure 12. Modal Amplitude F. vs time t for Traveling Haves - Unstable Case 
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Table 1. Comparison of Results for Showing Effects of Gas Dynamic 
Index (i) - (F x = 0, F x ' * 1, F 2 * 0, F 2 ' = 0, G x = 0, G^ = 0, G 2 = 0, 
G2* = 0) - Stable Cases (n - 60) - Standing Waves 



i 

X 

tl II 

i = 0 
K * 1 

t 

Exact 

Perturbation 

Exact 

Perturbation 


Solution 

Solution 

Solution 

Solution 

0.2 

0.19699 

0.18712 

0.19699 

0.18702 

0.4 

0.38335 

0.36426 

0.38336 

0.36386 

0.6 

0.55252 

0.52540 

0.55259 

0.52462 

0.8 

0.69856 

0.66518 

0.69885 

0.66400 

1.0 

0.81627 

0.77905 

0.81719 

0.77758 

1.2 

0.90132 

0.86340 

0.90354 

0.86186 

1.4 

0.95043 

0.91572 

0.95489 

0.91443 

1.6 

0.96159 

0.93461 

0.96936 

0.93399 

1.8 

0.93432 

0.91986 

0.94632 

0.92040 

2.0 

0.86985 

0.87244 

0.88656 

0.87466 

2.2 

0.77125 

0.79443 

0.79242 

0.79884 

2.4 

0.64330 

0.68895 

0.66783 

0.69602 

2.6 

0.49211 

0.56003 

0.51822 

0.57016 

2.8 

0.32469 

0.41247 

0.35021 

0.42593 

3.0 

0.14827 

0.25167 

0.17115 

0.26856 

3.2 

-0.03017 

0.08340 

-0.01142 

0.10366 

3.4 

-0.20430 

-0.08633 

-0.19032 

-0.06300 

3.6 

-0.36859 

-0.25158 

-0.35912 

-0.22566 

3.8 

-0.51829 

-0.40659 

-0.51242 

-0.37879 

4.0 

-0.64917 

-0.54604 

-0.64583 

-0.51727 

4.2 

-0.75738 

-0.66519 

-0.75583 

-0.63655 

4.4 

-0.83921 

-0.76006 

-0.83953 

-0.73278 

4.6 

-0.89118 

-0.82756 

-0.89450 

-0.80297 

4.8 

-0.91029 

-0.86556 

-0.91868 

-0.84504 

5.0 

-0.89444 

-0.87297 

-0.91049 

-0.85789 

5.2 

-0.84300 

-0.84979 

-0.86909 

-0.84143 

5.4 

-0.75726 

-0.79704 

-0.79483 

-0.79656 

5.6 

-0.64071 

-0.71679 

-0.68959 

-0.72514 

5.8 

-0.49885 

-0.61200 

-0.55708 

-0.62988 

6.0 

-0.33868 

-0.48646 

-0.40279 

-0.51427 

6.2 

-0.16794 

-0.34466 

-0.23365 

-0.38244 

6.4 

0.00574 

-0.19159 

-0.05742 

-0.23901 

6.6 

0.17554 

-0.03259 

0.11809 

-0.08892 

6.8 

0.33582 

0.12685 

0.28578 

0.06271 

7.0 

0.48221 

0.28126 

0.43974 

0.21080 

7.2 

0.61134 

0.42538 

0.57552 

0.35042 

7.4 

0.72042 

0.55435 

0.68999 

0.47703 

7.6 

0.80682 

0.66387 

0.78099 

0.58653 

7.8 

0.86778 

0.75032 

0.84693 

0.67550 

8.0 

0.90042 

0.81090 

0.88644 

0.74121 
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Table 2. Comparison of Results for F 2 Showing Effects of Gas Dynamic 
Index (i) - (F x a 0, F x * a 1, F 2 = 0, F 2 ' a 0, G x = 0, G^ = 0, G 2 = 0, 
G 2 ' » 0) - Stable Cases (n « 60) - Standing Waves 



i a 1 
K a 1 

i = 0 

K a 1 

t 

Exact 

Perturbation 

Exact 

Perturbation 


Solution 

Solution 

Solution 

Solution 

0.2 

0.00012 

-0.00485 

0.00003 

-0.00253 

0.4 

0.00113 

-0.01309 

0.00043 

-0.00938 

0.6 

0.00422 

-0.02223 

0.00205 

-0.01865 

0.8 

0.01060 

-0.02944 

0.00602 

-0.02784 

1.0 

0.02110 

-0.03215 

0.01336 

-0.03424 

1.2 

0.03582 

-0.02854 

0.02471 

-0.03553 

1.4 

0.05397 

-0.01795 

0.03497 

-0.03023 

1.6 

0.07375 

-0.00104 

0.05816 

-0.01798 

1.8 

0.09260 

0.02022 

0.07742 

0.00026 

2.0 

0.10749 

0.04283 

0.09521 

0.02232 

2.2 

0.11543 

0.06317 

0.10863 

0.04514 

2.4 

0.11407 

0.07764 

0.11493 

0.06516 

2.6 

0.10215 

0.08319 

0.11198 

0.07891 

2.8 

0.07988 

0.07788 

0.09878 

0.08357 

3.0 

0.04909 

0.06128 

0.07573 

0.07745 

3.2 

0.01309 

0.03461 

0.04478 

0.06031 

3.4 

-0.02375 

0.00071 

0.00929 

0.03352 

3.6 

-0.05653 

-0.03632 

-0.02639 

-0.00006 

3.8 

-0.08051 

-0.07164 

-0.05749 

-0.03640 

4.0 

-0.09180 

-0.10031 

-0.07945 

-0.07083 

4.2 

-0.08798 

-0.11801 

-0.08865 

-0.09863 

4.4 

-0.06850 

-0.12165 

-0.08296 

-0.11571 

4.6 

-0.03490 

-0.10993 

-0.06211 

-0.11921 

4.8 

0.00924 

-0.08353 

-0.02786 

-0.10793 

5.0 

0.05868 

-0.04515 

0.01610 

-0.08257 

5.2 

0.10711 

0.00077 

0.06458 

-0.04571 

5.4 

0.14796 

0.04864 

0.11144 

-0.00152 

5.6 

0.17534 

0.09236 

0.15041 

0.04468 

5.8 

0.18493 

0.12615 

0.17593 

0.08715 

6.0 

0.17466 

0.14534 

0.18393 

0.12042 

6.2 

0.14513 

0.14699 

0.17255 

0.14004 

6.4 

0.09957 

0.13035 

0.14243 

0.14314 

6.6 

0.04352 

0.09700 

0.09680 

0.12888 

6.8 

-0.01592 

0.05073 

0.0411 

0.09857 

7.0 

-0.07104 

-0.00295 

-0.01771 

0.05557 

7.2 

-0.11448 

-0.05745 

-0.07205 

0.00488 

7.4 

-0.14027 

-0.10594 

-0.11471 

-0.04743 

7.6 

-0.14458 

-0.14223 

-0.13988 

-0.09498 

7.8 

-0.12628 

-0.16158 

-0.14384 

-0.13189 

8.0 

-0.08716 

-0.16128 

-0.12556 

-0.15352 
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solutions are not exactly the same, the order of magnitude and behavior 
of results is similar, from Table 2, the same observations can be made 
for the behavior of F 2> There is, however, more error, quantitatively, 
between the results for exact and perturbation methods and a region of 
qualitative inaccuracy between the exact and perturbation solutions exists 
near t = 0. This takes the form of a difference in sign of F 2 between 
results from the exact solution as compared to the perturbations solution. 
This discrepancy occurred also in the other calculations performed (not 
shown) and will be discussed in more detail later in this chapter. 

hi Tables 3 and 4, a comparison of results is presented for modal 
amplitudes F^ and F 2 for a stable standing wave case. The initial 
conditions for the results in these tables are F^(0) = 0, F^'(0) = 1, 

F 2 (0) = 0, F 2 '(0) = 0, G 1 (0) = 0, G 1 ’(0) = 0, G 2 (0) = 0, G 2 '(0) = 0, 
n = 40, e = 0.1, and u = 0.1. However, these tables quantitatively 
present the effect of deviations of the ratio of the second acoustic 
frequency to the first from the integer value of 2 (this is controlled 
by the parameter K) . These results show that solutions for finite values 
of K are qualitatively similar to those for K = 0. This indicates that 
the ratio of the second acoustic frequency to the first does not have 
to be an integer in order to produce the type of behavior observed here. 

A ratio near an integer value will lead to similar results. Tables 3 
and 4 also allow a comparison to the results generated by the program 
in Appendix D for the approximate analytical solution (4.31). These 
results presented in the last column of Tables 3 and 4 can be compared 
to the fourth column in each of these tables to determine the accuracy 
of (4.31). These comparisons present further evidence that the neglect 
of gas dynamic nonlinearities does not have an important qualitative 


effect. 
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Table 3. Comparison of Results for Showing Effects of the Correction 
Variable (K) - (Fj a o, Fj' * 1, F 2 * 0» F 2 ' * 0 , Gj a o, Gj* ■ 0, 

Gg * 0, G 2 * 8 0) - Stable Case (n a 40) - Standing Waves 



i 

K * 

1 

1 

i a 
K a 

1 

0 

i = 0 
K = 0 

t 

Exact 

Solution 

Perturbation 

Solution 

Exact 

Solution 

Perturbation 

Solution 

Analytic 

Solution 

0.2 

0.19699 

0.18707 

0.19670 

0.19678 

0.19671 

0.4 

0.38335 

0.36396 

0.38172 

0.38210 

0.38186 

0.6 

0.55254 

0.52460 

0.54795 

0.54890 

0.54843 

0.8 

0.69867 

0.66358 

0.68903 

0.69093 

0.69029 

1.0 

0.81667 

0.77635 

0.79957 

0.60302 

0.80234 

1.2 

0.90247 

0.85935 

0.87537 

0.88124 

0.88075 

1.4 

0.95312 

0.91014 

0.91361 

0.92306 

0.92304 

1.6 

0.96699 

0.92744 

0.91310 

0.92740 

0.92820 

1.8 

0.94390 

0.91120 

0.87443 

0.89**68 

0.89666 

2.0 

0.88523 

0.86255 

0.80001 

0.82676 

0.83027 

2.2 

0.79388 

0.78374 

0.69393 

0.72690 

0.73221 

2.4 

0.67411 

0.67805 

0.56163 

0.59954 

0.60682 

2.6 

0.53129 

0.54969 

0.40948 

0.45016 

0.45947 

2.8 

0.37154 

0.40358 

0.24432 

0.28905 

0.29627 

3.0 

0.20133 

0.24520 

0.07300 

0.11102 

0.12386 

3.2 

0.02712 

0.08039 

-0.09782 

-0.06484 

-0.05085 

3.4 

-0.14491 

-0.08439 

-0.26190 

-0.23547 

-0.22097 

3.6 

-0.30901 

-0.24474 

-0.41335 

-0.39410 

-0.37989 

3.8 

-0.45992 

-0.39354 

-0.54662 

-0.53453 

-0.52152 

4.0 

-0.59287 

-0.52617 

-0.65660 

-0.65139 

-0.64055 

4.2 

-0.70357 

-0.63813 

-0.73872 

-0.74028 

-0.73261 

4.4 

-0.78821 

-0.72574 

-0.78928 

-0.79800 

-0.79446 

4.6 

-0.84362 

-0.78623 

-0.80583 

-0.82264 

-0.82407 

4.8 

-0.86746 

-0.81785 

-0.78760 

-0.81364 

-0.82069 

5.0 

-0.85849 

-0.81988 

-0.73571 

-0.77179 

-0.78490 

5.2 

-0.81682 

-0.79266 

-0.65317 

-0.69920 

-0.71852 

5.4 

-0.74409 

-0.73757 

-0.54464 

-0.59917 

-0.62456 

5.6 

-0.64352 

-0.65697 

-0.41585 

-0.47611 

-0.50704 

5.8 

-0.51967 

-0.55405 

-0.27311 

-0.33525 

-0.37090 

6.0 

-0.37816 

-0.43279 

-0.12272 

-0.18253 

-0.22172 

6.2 

-0.22516 

-0.29772 

0.02937 

-0.02427 

-0.06554 

6.4 

-0.06694 

-0.15382 

0.17767 

0.13305 

0.09140 

6.6 

0.09056 

-0.0063 

0.31714 

0.28307 

0.24291 

6.8 

0.24192 

0.13957 

0.44298 

0.41977 

0.38308 

7.0 

0.38242 

0.27867 

0.55061 

0.53776 

0.50650 

7.2 

0.50795 

0.40617 

0.63563 

0.63247 

0.60850 

7.4 

0.61493 

0.51774 

0.69408 

0.70029 

0.68527 

7.6 

0.70017 

0.60967 

0.72786 

0.73878 

0.73407 

7.8 

0.76090 

0.67898 

0.72020 

0.74667 

0.75327 

8.0 

0.79476 

0.72355 

0.68609 

0.72400 

0.74234 
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Table 4, Comparison of Results for F 2 Showing Effects of the Correction 
Variable (K) - (F 1 ■ 0, F A ' * 1, Fj • 0, F 2 ' ■ 0, G 1 • 0, G 1 * * 0, 
G2 s 0,G2' =0)- Stable Case (n * 40) - Standing Waves 



i * 
K = 

1 

1 

i * 1 

K = 0 

i = 0 
K = 0 

t 

Exact 

Perturbation 

Exact 

Perturbation 

Analytic 


Solution 

Solution 

Solution 

Solution 

Solution 

0.2 

0.00011 

-0.00400 

0.00015 

-0.00419 

-0.00192 

0.4 

0.00099 

-0.00996 

0.00136 

-0.01035 

-0.00696 

0.6 

0.00354 

-0.01599 

0.00479 

-0.01597 

-0.01337 

0.8 

0.00860 

-0.02013 

0.01141 

-0.01856 

-0.01885 

1.0 

0.01665 

-0.02071 

0.02149 

-0.01626 

-0.02113 

1.2 

0.02761 

-0.01669 

0.03443 

-0.00838 

-0.01856 

1.4 

0.04072 

-0.00793 

0.04866 

0.00438 

-0.01059 

1.6 

0.05458 

0.00479 

0.06180 

0.01995 

0.00208 

1.8 

0.06726 

0.01983 

0.07113 

0.03529 

0.01749 

2.0 

0.07663 

0.03491 

0.07416 

0.04696 

0.03278 

2.2 

0.08072 

0.04750 

0.06910 

0.05185 

0.04472 

2.4 

0.07803 

0.05521 

0.05549 

0.04793 

0.05038 

2.6 

0.06793 

0.05622 

0.03439 

0.03468 

0.04776 

2.8 

0.05081 

0.04957 

0.00841 

0.01343 

0.03626 

3.0 

0.02816 

0.03542 

-0.01367 

-0.01282 

0.01697 

3.2 

0.00244 

0.01507 

-0.04246 

-0.03978 

-0.00745 

3.4 

-0.02323 

-0.00917 

-0.05871 

-0.06265 

-0.03313 

3.6 

-0.04544 

-0.03425 

-0.06410 

-0.07695 

-0.05561 

3.8 

-0.06099 

-0.05679 

-0.05687 

-0.07938 

-0.07067 

4.0 

-0.06735 

-0.07356 

-0.03727 

-0.06854 

-0.07507 

4.2 

-0.06307 

-0.08196 

-0.00764 

-0.04530 

-0.06722 

4.4 

-0.04802 

-0.08040 

0.02778 

-0.01278 

-0.04759 

4.6 

-0.02355 

-0.06854 

0.06351 

0.02408 

-0.01872 

4.8 

0.00767 

-0.04745 

0.09363 

0.05926 

0.01509 

5.0 

0.04187 

-0.01945 

0.11279 

0.08668 

0.04845 

5.2 

0.07470 

0.01210 

0.11719 

0.10127 

0.07573 

5.4 

0.10173 

0.04330 

0.10529 

0.09987 

0.09206 

5.6 

0.11915 

0.07013 

0.07824 

0.08191 

0.09414 

5.8 

0.12425 

0.08905 

0.03976 

0.04959 

0.08095 

6.0 

0.11586 

0.09745 

-0.00434 

0.00764 

0.05395 

6.2 

0.09462 

0.09402 

-0.04711 

-0.03736 

0.01700 

6.4 

0.06286 

0.07894 

-0.08154 

-0.07807 

-0.02424 

6.6 

0.02438 

0.05389 

-0.10178 

-0.10757 

-0.06318 

6.8 

-0.01630 

0.02182 

-0.10409 

-0.12058 

-0.09336 

7.0 

-0.05322 

-0.01335 

-0.08754 

-0.11433 

-0.10952 

7.2 

-0.08237 

-0.04728 

-0.05423 

-0.08921 

-0.10854 

7.4 

-0.09956 

-0.07576 

-0.00908 

-0.04879 

-0.09001 

7.6 

-0.10237 

-0.09523 

0.04092 

0.00093 

-0.05639 

7.8 

-0.09018 

-0.10322 

0.08777 

0.05190 

-0.01267 

8.0 

-0.06428 

-0.09872 

0.12375 

0.09586 

0.03434 
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In Tables S and 6, a comparison of results are presented for modal 
amplitudes Fj and Fg for an unstable standing wave showing the effect of 
neglecting gas-dynamic nonlinearities. It can be seen that the gas 
dynamic nonlinearities have little qualitative effect on the results. 

In Tables 7 and 8, a comparison of results are presented for modal 
amplitudes F^ and Fj for an unstable standing wave case showing the effects 
of K. The results for zero and non-zero are qualitatively similar. 

These tables are representative of the cases that were investigated 
in the course of this research. Only cases involving standing waves were 
presented. The same behavior, however, can be observed for the cases 
involving traveling waves. 

In Table 9, a comparison of stability boundaries is presented 
based upon the interaction index (n) which is a measure of the strength 
of the combustion process. For standing waves and the given conditions 
shown, the stability limit for a process with gas dynamic nonlinearities 
considered and K = 0 is between 45-50. When both gas dynamic non- 
linearities and the correction variable are considered, the stability 
limit is increased to 67.5-69. Finally, when considering only the 
correction variable with no gas-dynamic non-linearity effect, the stability 
limit is 72-72.5. The results show that the neglect of gas dynamic 
non linearities slightly underestimates the stability boundary and that 
the increasing K increases the stability limit. 

In Table 10, a comparison of stability boundaries is presented 
based upon the interaction index for traveling waves. These results provide 
additional confirmation of the conclusions discussed in the previous 
paragraph and also illustrate the fact that standing waves are roughly 
twice as stable as traveling waves. This is consistent with the 


Table 5. Comparison of Results for Showing Effect of the Gas Dynamic 
Index (i) - (Fj = 0, Fj' * 1, F 2 * 0* F 2 * ■ 0, Gj * 0. G^' * 0, Gj * 0, 
G 2 ' = 0) - Unstable Case (n = 75) - Standing Waves 



i = 1 
K = 1 

i * 0 
K = 1 

t 

Exact 

Perturbation 

Exact 

Perturbation 


Solution 

Solution 

Solution 

Solution 

0.2 

0.19699 

0.18717 

0.19699 

0.18703 

0.4 

0.38335 

0.38215 

0.38336 

0.36402 

0.6 

0.55250 

0.52613 

0.55258 

0.52513 

0.8 

0.69847 

0.69244 

0.69881 

0.66518 

1.0 

0.81592 

0.80625 

0.81700 

0.77974 

1.2 

0.90030 

0.86736 

0.90294 

0.86535 

1.4 

0.94803 

0.92126 

0.95337 

0.91952 

1.6 

0.95671 

0.94184 

0.96695 

0.94089 

1.8 

0.92557 

0.92873 

0.94003 

0.92921 

2.0 

0.85572 

0.88276 

0.87584 

0.88529 

2.2 

0.75039 

0.80579 

0.77580 

0.81106 

2.4 

0.61484 

0.70079 

0.64411 

0.70941 

2.6 

0.45593 

0.57161 

0.48677 

0.58410 

2.8 

0.28147 

0.42287 

0.31118 

0.43963 

3.0 

0.09945 

0.25987 

0.12554 

0.28107 

3.2 

-0.08272 

0.08829 

-0.06194 

0.11388 

3.4 

-0.25864 

-0.08588 

-0.24371 

-0.05620 

3.6 

-0.42307 

-0.25669 

-0.41342 

-0.22350 

3.8 

-0.57177 

-0.41827 

-0.56608 

-0.38241 

4.0 

-0.70107 

-0.56512 

-0.69790 

-0.52773 

4.2 

-0.80745 

-0.69228 

-0.80594 

-0.63473 

4.4 

-0.88721 

-0.79547 

-0.88767 

-0.75934 

4.6 

-0.93641 

-0.87123 

-0.94056 

-0.83825 

4.8 

-0.95124 

-0.91703 

-0.96203 

-0.88902 

5.0 

-0.92862 

-0.93138 

-0.94961 

-0.91020 

5.2 

-0.86710 

-0.91381 

-0.90146 

-0.90121 

5.4 

-0.76771 

-0.86492 

-0.81711 

-0.86255 

5.6 

-0.63434 

-0.78635 

-0.69815 

-0.79563 

5.8 

-0.47368 

-0.68071 

-0.54872 

-0.70274 

6.0 

-0.29439 

-0.55149 

-0.37551 

-0.58703 

6.2 

-0.10598 

-0.40299 

-0.18717 

-0.45231 

6.4 

0.08258 

-0.22558 

0.00672 

-0.30299 

6.6 

0.26374 

-0.06815 

0.19691 

-0.14391 

6.8 

0.43198 

0.10715 

0.37555 

0.01979 

7.0 

0.58355 

0.28009 

0.53682 

0.18290 

7.2 

0.71613 

0.44478 

0.67707 

0.34025 

7.4 

0.82779 

0.59627 

0.79428 

0.48686 

7.6 

0.91635 

0.72909 

0.88737 

0.61813 

7.8 

0.97878 

0.83897 

0.95535 

0.72996 

8.0 

1.01111 

0.92224 

0.99664 

0.81891 
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Table 6. Comparison of Raaults for F2 Showing Effect of the Gas 
Dynamic Index (i) - (Fj * 0, F^' * 1, Fj ■ 0, Fj' = 0, Gj ■ 0, 

Gj/ s o t &2 s 0, 62* * 0) - Unstable Case - (n * 75) - Standing Waves 



i>l 
K » 1 

i ■ 0 
K * 1 

t 

Exact 

Perturbation 

Exact 

Perturbation 


Solution 

Solution 

Solution 

Solution 


0.00013 

-0.00548 

0.00003 

-0.00316 


0.00124 

-0.01544 

0.00054 

-0.01172 


0.00473 

-0.02691 

0.00257 

-0.02333 

0.8 

0.01210 

-0.03643 

0.00752 

-0.03484 

1.0 

0.02443 

-0.04077 

0.01670 

-0.04288 

1.2 

0.04197 

-0.03749 

0.03089 

-0.04728 

1.4 

0.06387 

-0.025518 

0.04994 

-0.03793 

1.6 

0.08806 

-0.005425 

0.07263 

-0.02258 

1.8 

0.11143 

0.02061 

0.0966 

0.00033 

2.0 

0.13027 

0.049007 

0.1186 

0.02813 

2.2 

0.14082 

0.07537 

0.13504 

0.05696 

2.4 

0.14005 

0.09513 

0.14236 ' 

0.08237 

2.6 

0.12628 

0.10426 

0.13796 

0.09995 

2.8 

0.09967 

0.10004 

0.12064 

0.10606 

3.0 

0.06242 

0.08152 

0.09109 

0.09851 

3.2 

0.01860 

0.04984 

0.05192 

0.07688 

3.4 

-0.02636 

0.008166 

0.0076 

0.04283 

3.6 

-0.06636 

-0.03858 

-0.03628 

-0.00096 

3.8 

-0.09548 

-0.08436 

-0.07362 

-0.046816 

4.0 

-0.10889 

-0.12285 

-0.09877 

-0.09134 

4.2 

-0.10354 

-0.14826 

-0.10733 

-0.12758 

4.4 

-0.07873 

-0.15626 

-0.09688 

-0.15016 

4.6 

-0.03636 

-0.14461 

-0.06742 

-0.15524 

4.8 

0.01909 

-0.11361 

-0.02153 

-0.14104 

5.0 

0.08104 

-0.06617 

0.03583 

-0.10827 

5.2 

0.14148 

-0.00758 

0.09778 

-0.06010 

5.4 

0.19206 

0.055118 

0.15638 

-0.00186 

5.0 

0.22528 

0.11399 

0.20358 

0.05954 

5.8 

0.23565 

0.16121 

0.23235 

0.11651 

6.0 

0.22067 

0.19016 

0.23779 

0.16168 

6.2 

0.18130 

0.19633 

0.21796 

0.18888 

6.4 

0.12201 

0.17793 

0.17434 

0.19397 

6.6 

0.05019 

0.13633 

0.11182 

0.17547 

6.8 

-0.02482 

0.07595 

0.03808 

0.33477 

7.0 

-0.09302 

-0.02526 

-0.03737 

0.07614 

7.2 

-0.14503 

-0.07144 

-0.10448 

0.006209 

7.4 

-0.17338 

-0.140235 

-0.15407 

-0.06680 

7.6 

-0.17347 

-0.19368 

-0.17900 

-0.13401 

7.8 

-0.14428 

-0.22450 

-0.17519 

-0.18700 

8.0 

-0.08859 

-0.22808 

-0.14216 

-0.21888 
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Table 7. Comparison of Results for Showing the Effects of the 

Correction Variable (K) - (Fj * 0, F^ * * 1, Fj « 0, F 2 * * 0, Gj * 0, 
Gj 1 s 0, Gj * 0, G 2 ' * 0) - Unstable Cases (n * 70) • Standing Waves 




1*1 
K * 1 

i * 
K * 

1 

0 


i * 0 
K * 0 

Exact 

Perturbation 

Exact 


Perturbation 

Analytic 

Solution 

Solution 

Solution 


Solution 

Solution 

0.19699 

0.187155 

0.19670 


0.19687 

0.19675 

0.38335 

0.36462 

0.38172 


0.38261 

0.38217 

0.55251 

0.52587 

0.54791 


0.55028 

0.54942 

0.69850 

0.66615 

0.68878 


0.69371 

0.69249 

0.81604 

0.78072 

0.79865 


0.80767 

0.80630 

0.90066 

0.86596 

0.87280 


0.88811 

0.88696 

0.94887 

0.91928 

0.90775 


0.93226 

0.93184 

0.95842 

0.93925 

0.90166 


0.93878 

0.93468 

0.92863 

0.92554 

0.85477 


0.90772 

0.91061 

0.86067 

0.87904 

0.76962 


0.84061 

0.84612 

0.75771 

0.80166 

0.65109 


0.74802 

0.74903 

0.62481 

0.69648 

0.50590 


0.61107 

0.62333 

0.46860 

0.56735 

0.34194 


0.45807 

0.47407 

0.29658 

0.41903 

0.16733 


0.26849 

0.30714 

0.11649 

0.25679 

-0.01027 


0.10616 

0.12906 

-0.06441 

0.08639 

-0.18416 


-0.05313 

-0.05327 

-0.23973 

-0.08618 

-0.34860 


-0.25962 

-0.23280 

-0.40414 

-0.25437 

-0.49854 


-0.42949 

-0.40262 

-0.55328 

-0.41415 

-0.62912 


-0.58151 

-0.55621 

-0.68301 

-0.55829 

-0.73533 


-0.70960 

-0.68766 

-0.78998 

-0.68250 

-0.81197 


-0.80859 

-0.79192 

-0.87041 

-0.78258 

-0.85399 


-0.87449 

-0.86494 

-0.92055 

-0.85525 

-0.85736 


-0.90458 

-0.90388 

-0.93689 

-0.89811 

-0.82008 


-0.89755 

-0.90718 

-0.91671 

-0.90981 

-0.74306 


-0.85356 

-0.87463 

-0.85886 

-0.89005 

-0.63045 


-0.80213 

-0.80736 

-0.76447 

-0.83962 

-0.48925 


-0.66264 

-0.70785 

-0.63724 

-0.76029 

-0.32818 


-0.52309 

-0.57982 

-0.48340 

-0.65486 

-0.15634 


-0.36111 

-0.42807 

-0.31100 

-0.52691 

0.01809 


-0.18308 

-0.25635 

-0.12888 

-0.38083 

0.18852 


-0.07645 

-0.07712 

0.05445 

-0.22156 

0.34995 


0.10764 

0.10864 

0.23170 

-0.05445 

0.49844 


0.37422 

0.29176 

0.39723 

0.11485 

0.63023 


0.46032 

0.46509 

0.54707 

0.28072 

0. 74089 


0.69003 

0.62177 

0.67846 

0.43759 

0.82500 


0.81052 

0.75551 

0.78917 

0.58027 

0.87632 


0.89873 

0.86085 

0.87687 

0.70394 

0.88888 


0.95061 

0.93330 

0.93868 

0.80457 

0.85847 


0.96357 

0.96962 

0.97109 

0.87877 

0.78419 


0.93647 

0.96786 
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Table 8. Cougar is on of Rssults for f2 Showing the Effects of the 
Correction Variable (K) - (Fj » 0, F^' ■ 1, F 2 ■ 0, F 2 ' ■ 0, Cj ■ 0, 
G^' r. 0, Gg ■ 0, 62' * 0) - Unstable Cases (n ■ 70) - Standing Waves 



| * ■ 

1 

1 

i ■ 
K ■ 

1 

0 

i « 0 
K • 0 

t 

Exact 

Perturbation 

Exact 

Perturbation 

Analytic 


Solution 

Solution 

Solution 

Solution 

Solution 

0.2 

0.00013 

-0.00527 

0.00017 

-0.00562 

-0.00336 

0.4 

0.00120 

-0.01466 

0.00165 

-0.01557 

-0.01219 

0.6 

0.00456 

-0.02534 

0.00619 

-0.02603 

-0.02343 

0.8 

0.01160 

-0.03411 

0.01544 

-0.03276 

-0.03305 

1.0 

0.02332 

-0.03789 

0.03024 

-0.03218 

-0.03709 

1.2 

0.03992 

-0.03451 

0.05013 

-0.02232 

-0.03264 

1.4 

0.06057 

-0.02299 

0.07306 

-0.01862 

-0.01865 

1.6 

0.08330 

-0.00396 

0.09555 

0.021911 

0.00367 

1.8 

0.10517 

0.02047 

0.11324 

0.04915 

0.03093 

2.0 

0.12271 

0.04692 

0.12170 

0.07264 

0.05810 

2.2 

0.13242 

0.07126 

0.11751 

0.08675 

0.07947 

2.4 

0.13150 

0.08922 

0.09917 

•0.087046 

0.08978 

2.6 

0.11839 

0.097153 

0.06770 

0.071414 

0.08535 

2.8 

0.09328 

0.11651 

0.02679 

0.04063 

0.06502 

3.0 

0.05822 

0.07467 

-0.01761 

-0.001408 

0.03053 

3.2 

0.01701 

0.04468 

-0.05827 

-0.04828 

-0.01346 

3.4 

-0.02528 

0.00567 

-0.08789 

-0.09196 

-0.06005 

3.6 

-0.06295 

-0.03775 

-0.10041 

-0.12426 

-0.10122 

3.8 

-0.09049 

-0.07994 

-0.09221 

-0.13827 

-0.12919 

4.0 

-0.10337 

-0.11504 

-0.06290 

-0.12986 

-0.13785 

4.2 

-0.09872 

-0.13779 

-0.01566 

-0.09854 

-0.12403 

4.4 

-0.07590 

-0.1442 

0.04299 

-0.04785 

-0.08825 

4.6 

-0,03665 

-0.13261 

0.10404 

0.01493 

-0.03490 

4.8 

0.01487 

-0.1032 

0.15731 

0.07988 

-0.02828 

5.0 

0.07254 

-0.05894 

0.19315 

0.13602 

0.09131 

5.2 

0.12595 

-0.00478 

0.20424 

0.17311 

0.14357 

5.4 

0.17637 

0.05266 

0.18706 

0.183575 

0.17558 

5.6 

0.20783 

0.10617 

0.14284 

0.163922 

0.18070 

5.8 

0.21824 

0.14862 

0.07764 

0.11562 

0.15640 

6.0 

0.20519 

0.17409 

0.00147 

0.04509 

0.1049b 

6.2 

0.16949 

9.17865 

-0.0733 

-0.03703 

0.03331 

6.4 

0.11514 

0.16089 

-0.13397 

-0.11757 

-0.04784 

6.6 

0.04883 

0.12229 

-0.16972 

-0.182899 

-0.12564 

6.8 

-0.02090 

0.06705 

-0.17346 

-0.22122 

-0.18709 

7.0 

-0.08483 

0.00165 

-0.14308 

-0,22468 

-0.22125 

7.2 

-0.13429 

-0.06594 

-0.08205 

-0,19089 

-0.22111 

7.4 

-0.16228 

-0.12723 

0.00095 

-0.12367 

-0.18495 

7.6 

-0.16442 

-0.17433 

0.09312 

-0.03255 

-0.116J0 

7.8 

-0.13960 

-0.20087 

0.17949 

0.06839 

-0.02652 

8.0 

-0.09020 

-0.21857 

0.24522 

0.16276 

0.07753 
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Table 9. Comparison of Stability Bound arias Baaad on tha In t tract ion 
Indax (n) - « 0, Fj* * 1, F 2 • 0, F 2 ' ■ 0, 0! ■ 0, G*’ ■ 0. 

G 2 * 0, G 2 ' * 0) - Standing Havas - Epsilon - 0.1 


Gas Dynamic Indax 
Corrtction Variabla 


n 


i * 1 
K ■ 1 


i * 0 
K ■ 1 


Stability Boundarias 

Exact Solution Perturbation Solution 

- Stable - Unstable n - Stable - Unstable 


07.5 - 69 


67.5 - 69 


72 


72.5 


72.5 


73 


i * 1 
K « 0 


45 - 50 


45 - 50 
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Table 10. Comparison of Stability Boundaries Based on the Interaction 
Index (n) * (F^ ■ 0, Fj* ■ -1, F 2 * 0, F 2 * * 0, Qj * 1, Oj' ■ 0, 

62 1 0, Gj' ■ 0) - Traveling Haves - Epsilon - 0 .’ 
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approximate analytical stability equations (4.31) and (4.52). The 
perturbation method tends to predict slightly higher stability limits 
than the exact solution method for both standing and traveling waves. 
Within the accuracy of the tabulated values, this is apparent only in 
the first two rows of Table 10. 

In Table 11, a comparison of the effect of different initial 
conditions imposed on the stability boundaries for both standing and 
traveling waves is presented. From the results of two sets of initial 
conditions for each case, it can be seen that the varying of initial 
conditions has no significant effect on the stability boundaries for 
both standing waves or traveling waves. 

In Table 12, the variation of the stability limit with e is 
presented for standing waves. From Table 12, the results show that the 
smaller the term epsilon the greater the stability limit. Therefore, 
the order term has a significant effect on the interaction index. In 
Chapter 4, a relation was proposed for the case of i = 0 and K = 0 
which was n = C/e where C is a constant. Assuming the validity of the 
relation, the values for this constant are given for each given epsilon 
and interaction index. This shows that, in general, C is a weak function 
of e. 

In Table 13, a comparison of the effect of e is presented for 
traveling waves when both gas dynamic nonlinearities and correction 
variables are considered. Again, the results show that the smaller the 
term epsilon, the greater the stability limit. The perturbation method 
again predicts slightly greater stability limits than does the exact 
solution method. Therefore, again, the order term has a strong effect 
concerning the stability of combustion. 


Table 11. Comparison of the Effect of Different Initial Conditions 
Imposed for Standing and Traveling Waves for i * 1 and K = 1 
Epsilon * 0.1 

(a) Standing Waves - 1. Fj * 0, F^' = 1, F 2 * 0, F 2 ' s 0 

G i =0, Gj* = 0,G 2 =0, G 2 B 0 


Initial Condition 
Sets 

1. 


2 . 


2. F A = 1, Fj = 0, F 2 = 0, F 2 = 0 

Gj = 0 , G^ ' = 0 , G 2 = 0 , G 2 * = 0 

Stability Boundaries 


Exact Solution 
n - Stable - Unstable 

67.5 - 69 
65 - 70 


Perturbation Solution 
n - Stable - Unstable 

67.5 - 69 
65 - 70 


(b) Traveling Waves - 1. 


*1 

Gl 


0, F ! * = -1, F 2 = 0, F 2 ' = 0 


1 1 * = 0 , G 2 = 0 1 G 2 = 0 


i _ 


= 1, F 1 ' = 0, F 2 = 0, F 2 ' 


G^ = 0 , G ^ ' - -1, G 2 = 0, G 2 


= 0 


Initial Condition 
Sets 


Stability Boundaries 


Exact Solution 
n - Stable - Unstable 

27.5 - 28 


Perturbation Solution 
n - Stable - Unstable 

31.5 - 32 


2 


27.5 - 28.5 


31 - 31.5 
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Table 12. Comparison of the Effects of the Order Term Epsilon 

(Fj * 0, Fj = 1, F 2 * 0j F 2 - 0| Gj 3 0, Gj 1 * 0, G 2 3 0* 

G'i' 3 0) - Standing Waves - when i * 1, K =1 


Epsilon 

Stability Boundaries 

Exact Solution 
n - Stable - Unstable 

Perturbation Solution 
n - Stable - Unstable 

0.05 

107.5 - 110 

107.5 - 110 

0.1 

67.5 - 69 

67.5 - 69 

0.2 

48.5 - 49.5 

48.5 - 49.5 


Constant 
C = nc 














Table 13. Comparison of the Effects of the Order Term Epsilon (F 
F 2 = 0, Gj = 1 , G. ' = 0, G 2 = 0, G_* = 0) Traveling Waves - when 
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Thus, from these representative tables of results , it is observed 
that the correction variable is important in the stability of standing 
waves, but does not play a major role in the stability of traveling waves. 
It is observed that the gas dynamic nonlinearities seem to have little 
influence on the stability of either standing or traveling waves. It 
is observed that initial conditions of the modal amplitudes have little 
or no influence in the stability of either standing or traveling waves. 

And finally, it is observed that the order term epsilon and, the inter- 
action index governing the strength of combustion in the process are 
strongly coupled thus affecting the limits of stability. 

Before completing this chapter, it is desired to investigate the 
sign discrepency mentioned previously between the exact and perturbation 
solutions for f which occur near t = 0. For simplicity, it will be 
assumed that i = K = 0 and that for t « 1 the first modal amplitude can 
be represented with sufficient accuracy by f^ = sint. Then, the 
equation for f 2 will be solved and the result simplified for t << 1. 

This will be done first for u = 0 and then for u i 0. For u> = 0, 

(3.21) leads to 


d2f r 1 

+ 4f 2 = h eun 1 1 - cos2t 


(5.1) 


with initial conditions 


f 2 (0) = 0 


f 2 (0; = 0. 
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Evaluating the homogeneous and particular solutions by the usual manner 
and evaluating the constants, the results become 

f 2 = Yg cun £l - cos2t - t sin2tj . (5.2) 

In terms of the perturbation parameters (4.1), equation (5.2) can be 
written as 



am e(l - cos2£) - n sin2£ 


» • 


(5.3) 


To the order of approximation e which the perturbation solution should 
model, equation (5.3) becomes 


f 


2 


1 _ 

16 


umn sin2C + 0 (e). 


(5.4) 


By expanding equation (5.2) into a Taylor series expansion of three terms, 
equation (5.2) becomes 


1 - 4 

f 2 = 24 eunt + • • • (5.5) 

which is always positive. 

Therefore, the exact method for small time will yield f 2 modal 
amplitude always as a positive quantity. 

By imposing identical conditions to the perturbation equations 


(4.12), the result becomes 
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(5.6) 


with the condition 


b 2 (o) = 0. 

Solving equation (5.6), 

B 2 = - unn. (5.7) 


Recalling that f 2 = B 2 sin2£, the result becomes 

i _ 

f 2 = - wnn sin2£ + 0 (e) 


(5.8) 


which is identical to the result of equation (5.4) for the wave equation 
solution. Thus, the perturbation method gives the correct result. It 
can be seen that for t << 1 the exact solution predicts a positive f 2 
and by inspection of equation (5.8), the perturbation method predicts a 
negative f 2> This is precisely the behavior observed in the numerical 
solutions. 

For us ^ 0, a similar analysis can be performed. The appropriate 
equation for f 2 is now 


d 2 f df 
d^ + “ dt* 


+ 4f„ 


kieionj^ 


1 - cos2t 


(5.9) 


with conditions 


f 2 (0) ■ 0 
f 2 ’(0) * 0. 

Solving the homogeneous and particular solution by the usual manner and 
evaluating the appropriate constants the result becomes 




/IT 


- 52 - 


t + 


1_ 

16 


eun 


— en sin2t 

O 


=]. 


(5.10) 


Expanding (5.10) for small u into the appropriate Taylor series, expanding 
and neglecting terms of 0(u>) leads to 



eun 


[‘ 


cos2t - t sin2t 


(5.11) 


which is identical to (5.2). 

By imposing the identical conditions on the perturbation equation 
(4.12), the resulting equation become 



+ h oB 2 



(5.12) 


with the condition 
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b 2 (o) « 0. 

Solving equation (5.12) by the usual manner* evaluating the constants, 
and transforming the perturbation variables to real time variables 

#J . In St. (5.13) 

This is always negative for t « 1. Expanding the exponential function 
by the Taylor series expansion and neglect terms of o(w) leads to 

f 2 * n S * n (5.14) 

which is identical to (5,8). 

To observe the behavior of equation (5.10) for small time, 
expand this equation into a Taylor series of 0( t 1 * ) . Expanding and 
grouping terms according to their order of magnitude, the terms of 
0(1), 0(t), 0(t 2 ), 0(t 3 ) vanish. Therefore, fj is comprised of terms 
from 0(t 4 ) which is 


- _ eniSt 4 T. . 3uS 2 . w 4 1 

2 " 24 L 1 8 64 J . (5.15) 

4 

Again, for any small time t, f 2 is always positive since t is always 
positive. Neglecting higher powers of u, the resulting equation becomes 
equation (5.5) for the undamped case. Again it can be seen that the 
exact and perturbation methods predict opposite signs for when t << 1. 
These results are based on approximations and cannot be considered 


definitive. They do, however , lend plaueibility to the numerical results 
discussed earlier. It is believed that this sign discrepancy is due 
to the inability of the perturbation solution to accurately represent 
the exact solution for t « 1 and not due to any error .In the computer 
program used to compute the perturbation solution. 



Chapter 6 


CONCLUSION AND RECOMMENDATIONS 

The primary objective of this presentation has been the development 
of analytical techniques to solve the problem of combustion instabilities 
occurring in an annular combustion chamber . The analytical techniques 
used were the modified Galerkin method applied to the acoustic wave 
equations which yielded a set of time-dependent modal amplitude equations 
and the two-variable perturbation method which yield a set of time- 
dependent equations which approximated the behavior of the first set of 
equations. Both methods produced results which were relatively easy to 
apply and used the Runge-Kutta algorithm which required little computation 
time. An alternative approach to solve this problem would be a finite 
difference approach. However* difficulties can be foreseen in the 
development of the finite difference equations modelling the problem 
along with the complications occurring due to the boundary conditions of 
the problem. Thus, the benefits of the methods discussed in this thesis 
can be appreciated. 

From the numerical and graphical presentation of results in Chapter 
5, the following observations can be made. First, the effect of the gas- 
dynamic nonlinearities seems to be small in both methods of analysis for 
velocity sensitive combustion. This point can be observed from a 
quantitative comparison of the tabular results or by observing the effects 
of this condition on the stability boundaries. Second, the effect of the 
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correction variable modelling tha physical boundaries of the chamber seems 
to have a significant effect in both methods of analysis for velocity 
sensitive combustion. By including the effect of this correction variable* 
a significant increase occurs in the interaction index which is the 
criteria for the stability of the system. However, this effect seems to 
be more significant for the standing wave case than the traveling wave 
cases. The effects of initial conditions for the time dependent equations, 
the numerical value for the burning rate and step size of integration, 
seem to have very little significance in the measure of the stability 
limits of velocity sensiti e combustion. However, the order term epsilon 
has a strong effect upon .he stability of the problem. This is to be 
excepted since the order term is the measure of the effect of non- 
linearities occurring in the system. The increase in this value corresponds 
to a decrease in the stability limit which is physically reasonable. 

In this study, the effect of time delay of the combustion process 
was neglected. However, time delay has been found in other studies to 
be an Important phenomena in correctly modelling t».e actual problems of 
velocity sensitive combustion. It is recommended that this effect can 
be incorporated by including the corresponding terms with j * 1 in the 
acoustic wave equations (3.20). A corresponding set of perturbations can 
then be derived to account for time delay and both these equations and 
equations (3.20) can be numerically evaluated by modlfing the existing 
Runge-Kutta programs presented in the Appendices. It is also recommended 
that an experimental program be developed to measure the effects of 
velocity sensitive combustion in an annular combustion chamber. Once 
achieving this goal, one could correlate the measurement results to the 
analytical results that have been presented to ascertain the validity of 
this analysis. 


Since instability of combustion is sensitive to smell changes in 
engine geometry and operating conditions , a particular engine must be 
subjected to a large number of firings before its deciders can say 
confidently that it is free from instability. Vith a luge engine such 
testing can account for a substantial part of development costs. Herein 
lies the importance of devising reliable theories of instability and 
inexpensive tests of a propellant's acoustical characteristics. Until 
instability of combustion is understood well enough so that it can be 
eliminated while an engine is in the design stage, rocket engines must 
continue to be intensively tested for stability— particularly when 
the lives of astronauts will eventually depend on safe, reliable 
operation of the engine [17]. 
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t - 5 < T 

t -*T < € 

t - T > T . (A 

Graphically representing equation (A. 2) results in Figure A2. 

J(t - o| 


0 t - C t 

Figure A2. Step Time Delay Function J(t - 5) 

Substituting into the general time-delay integral the particular step 
function in terms of the non-dimensional variable 5 


or 


J(t - C) s 


t 


J(t - o = 


ft - T du ft dcu 

«ar* Mar 4 *- (A 

0 t-T 

Therefore, simplifying equation (A. 3) 

u(t) = u Q (t) - Wg ( t - t) (A 

where w Q (t) is a generalized function of time and u^ft - t) is 
functional time delay. 
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2/2. .32* 52/2. #fip 
2/2 ..Gj. j ?/?. ,fip 

2/2..32452/2».Fl° 

2/2 . » G 2 4 52/2, .rip 


1«TJ/’.*F2?*U2/2.»GIP< V2/2.#Gi!t>«W2/>. > 
M3M*£9(T4H/2«4Fi4 P2/2.*F 2* 82/ 2, # G1 4* 2/ 2 . *G2* 32/2 . #F l P 
l4T2/’.»F2»4U2/3.»GlP«V2/2.»G2P®M2/2.) 

P4»^*£l(t« H»Fi< P3#F2-* 3 3»Gl*S3#G2*S3#Flp*r 3»F2®*U3»GtP 
14 v3»32?4*3) 

C**H*c2(T4H#Fl4P3»F24as#GHRS^G2 4S3#Fl®4 T3»F2P4U5#G1P 
1 ♦ V3 »12P4 M3 ) 

«4*H*E3(T4H*Fl4P3»F2433#Gl4ft3»G24S3»FlP4T 3#F2?4U3/G1® 
14^3»G2P4H3 ) 

54«H*£<*(T« H#Fi4P3»F24 33,Gt4R3»G24S3#Fl»4 T3#F2P4U3»Gl» 
14w3*G2P4M3 > 

T«*H*£5(T4H,F14P3*F24 33#G14P3#G24S3» : ‘1P4 r 3#F2?4U3,G1P 
1 ♦ »3* 32P 4 <«3 ) 

U<.»H*£GCT4 H»fl *P3» F2433»G14P 3 »r.24S3» FI® 4 T3»F2P* U3* Gl» 
14V?#32°«^3) 

^;*H*£?CT4H»F1 4 p 3, F24(J3»314R 3#G 24 5 3 *F1P4 T3»F2?4U3#G1» 

1 4 43 » G2P« M3 ) 

M4*H*£e(T4 H»FI*P3#F24 3 3»G14P3#G24S J»"l®4 T 3»F2P4'j3»GiP 
1* »3/ 12®«*3> 

FI«Fl4(f»i #j.*P2»2 # #P34P4)/6. ) 

F2»* 24 ((31 4 2.0 24 2.*.3S«3 4 )/o. > 

51®G14C(PI«2.*P2 42.*R34K4 >/S. ) 

£ 2«3 24 C( Si 4?. * 5 24 2. •$!«:»*)/ 6. ) 

Fl a ®FlP4((n42.*T2*2.*T34T4)/S.) 

F2?«F’P*C(Ul4’.*U242.*U3«U4)/6.) 

Gif’* J1P4(( m«?.* /2*2. 4 Z3« V4 )/ S. ) 

«G2®4 ({ Wl*2.*W242.*i*34 M4)/i. > 

1C’! ®-Fl®4 £PS*C "Fl *F2 • G i *r,24 1.5*{ Gl?»J2P4 F1P*F2?) ) 
iC®2 = -F2®4 £®S*{ -0.25*( ji*Gl- r l*Fi)*0. :5®(FIP*F1P-GIP 
1 *Gi i* ) ) 

4 CP 3 ■-31®4 E®$*( -Fi*G’« F2*G140. 5*( F IP * 1 2® *G l? *F 2P> ) 
4C»»4«-a2*« t B S*f 0.5*F l*3i*0.5* r l°*Gl P) 

r*T4H 

L=L‘l 

if a .ei. 3) io ro no 

3G TO 100 

110 *R1T£ (3*50) T •Fi»F2#3i»32»ASPl»«CP2»4CR3»ACP4 
50 FCPHAT (lX*F3.b*3X»Fl j. 3# 2X»R1).S> 3X» F10.5*2X»F10.S»SX 
i»Fn.i*2X»Fl9. J» 3X»F10.5»2X#F10.5> 

L = l 

101 CCNTIMUc 
CALL EXIT 
tSO 



ooo n n o n n ri n n n n r. n o n^ooo 


GA«Y H MCOONAlf) NON-IINFA* C0M3UST ION £ 9U AT I ON'i 0ERI/S3 
M3H THE 

TWJ-¥Afl!A8L£ PERTURBATION HETHDO 

GOVERNING EjUATUNS-CASE - NO TIME OSLAY CJ»0) 

Elf T» 4l#3l»A2#32»A3»33#A4#94)*(-0.5*SlG3*Al-0.S*«*3l 
I-C. 5*1 •< 4 1 *4 3*9i*3 3<4?*A4*32*a4>-0.5*N*.<3AP*( 31*A3-Al 
2* 63*32*A4-A2*34 )>*£PS 

E 2 C T # Al# 31 » A2, 32»A3#33»A4#34) ■( -0.5*5 1G3* flU 0. 5** *A l 
t-3.5*I*(Ai*9S3l*A3-A4*324A2*J4)-J.5*N*43AR*<.\l*AS*31 
2*53«A2*A4«3?*34 ))*£P$ 

£2<T, A 1,31 , A2,*2#a 3#3 3# A 4*3 4) •< * 0. 5*$ 13 3 * A 2* 0. 5*< *3 2 
t* J.5*I*<Al *44« ? 1*94-4 2*4 S*32* J3)* w5*N*43AP*( 44*31-Al 
2*£4M2*33-A3*92 > )• £®S 

E4(T# Al#3l # 42# 3 2# 4 3#3 3# A 4, 34) *t*2. S*S 1S3*32« 0.5 *,<*a 2 
l-;.S*l*(14*4U4*3i4A3-32-*2*33)-)o*N*43AR*<4l*A4'9i 
2*c4"42*Ai*9?*33) )*£PS 

£!( T# Al#31 # A2# 3 2# A 3#43#44*34) *( • 0. 5*5 133* A 3*4. *<<*3 3 

1- S-125*! *( 4 2*4 2*J2*’J?«3l*31-U*4l,-0. 1>5*W34P*N*< Al*3I 

2- 42* 32 )>*£P$ 

£6 C T # 41*51 # 12# 32#A3#i3#A4*94) ■ ( - •). 5*S 13 3 *3 3« 4. ** *a 3 
I-C - 25*1*14 2*'92At*3t)-0.G625*4 3AR*N*(I2* A2-92*12-U*A1 
2 ♦ El * 3 1 ) )*£ PS 

£7<T# A|» 31 #4 2, ?2#A3#j3#A4*34) *( - 0. 5* $ 1 3 3 * A 4-4. *x *94 
X-C-25*! *(3 1*32» l*A2)-').12‘s*W34«*N*( Al *32* 42 *31 ) )*£®5 
£J(T*4l#3l#42»32#A3#3S#44#34)»(»0,0*5Jj3*9 4«4.»i(*44 
A* C. ?5 *1 *t J2*4l« A 2* 31 )«0.l2S*<SjAfi*,N*( Al* A 2*3 1 • 3 2 ) ) *£P $ 

DIMENSION Ft(5<)'5>»F2<5D<5)#Gl( 500>#32( 530) 

3 • m i. H # N • i 

L«i 

r»o. 

6£AJ INITIAL CONDITIONS *NO CONSTANTS 

SPS-OSOE* TERM (EPSILON) 

VMAfc-STEA'JY STATE $U«NING fi»T£ 

A* I NT £4 AC T ION JNOEX 

H-CORPECT ION VARIABLE ( 3 AFFLES# WALL LI NI N0$ # NOZZLES 
»£TC. ) 

1-045 DYNAMIC INDEX 

Al»31# A2* 3?#A3#J3# A4«i4*H0i)AL AMPLITUDES SDEFF I Cl ENT S 

READ ( 5 # 1 •) ) «»*J4F #N» K* NN Ax# EPS# I 
1C FQP«4! C Sri0.4#I j#2Fl).4) 

*£A 3 (SilM Al# )l»A2#32#Aj#3J» V4 * j 4 
12 PC-MAT C!»r 13. A ) 


P-A3 STEP SIZE 


REPM®UC1BIUTY OF THE 
ORIGINAL PAGL lb PGOR 
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*•43 (5,14) H 
U FC**N4T (F10.5) 

5m«N1AX/E*S 

nPITi <&#20) SIGd**3AUN#K#H#£PS*I 
2U FCS^4T <IX#4X#* SIGMA i4M#Fl0.4#/#l X#4X# » STEADY STATE 
1 SONINS 9 4T£-(HS4F> , #*‘l0.4#/#tx#4X#*lMTES4CTlQ« INDEX 
2-CN ?* #Fl9.4#/#lx#4X#* CO® SECTION VAX I I 3LE-< X) UFIO.A#/ 
3#lX»4X#*SrCP S)ZE-(H)'*Fi:.4#/#lX#4X#'£PSIL0N-(EPS>* 
4#F10.4#/#IX*4X#*3AS DYNAMIC INDEX- ( 1 )»# FI 0. 4# // ) 

AA1TC (6#3C) *1*31*42*32# AS#3 3*44*34 
30 F On MAT (IX#* INITIAL C3N0 1 TI ON S* #// # IX # 4X # • A l • #F 13 .4# A X 
*•*31* #F1Q. 4*/*lX#4X** 42 *» FI 0.4# 3x#*92'*ri0.4#/*lX*4x 
2* *4 S* #F 13. 4#3X# ' 33' *"10.4*/,IX#4<# « A4 »#Ft 0. 4, 3X#*34* 
3*Fi:.4#//) 

*"I T E (3*4 0) 

40 F C M 4 I ClX*TX»*TjMfc**?',:. *41*# )X#*&1 *#12X# • 4 2* # l IX# * 32 ’ 
l#i:x» • A$*. UX# *3J*#1DX# *4»**| Cx#*i}4*#/> 
c 

C 4 1N-|£ X JTTA ALC3* I Thm 
C 

DC 170 J * l#N*AX 
C 

Pl*3*ElU# 41*31# A2# 32# 43*3 3* 4 4, 34) 

JIM *E2(T. A l #31 #42# 32# 43* 33# 44*34 ) 

P1«H*£!(T# At #31*42*32 #41*33*44* J4) 
3l*H*£4(T#Al#31#A2*d2*43*3J*4<'*?4) 

Tl*H*£3(T#*l»Jl# 42*32# A3* 33*44*34) 

U I *H* £ 3( T# A l # j 1 # 4 2 # 3 2 # A 5 *3 3 * 4 4 # 34 ) 
YI«'3*ET(T#a;#31#A2#8!#A3*3J*44#34) 

#l*H*E3(T* 41*31# A 2# 32# 4 J # J 3*44*24) 

C 

P2«M*£l(T*H/2.»4i'?i/2, # )l4 si / ?.,A’«x 
t*Tl/2.#3S* , w ,l/2.#44*Vl/?.#J4**l/2.) 

: 2 *H*E 2 (T* H /2 .* 41 4 PI /! . * 31 « 31 />., i 2 *s 1 /> i?# >i /*.,* 3 

i« 1 1 / 2 • #3 3# Jl/2.#44« /* /?. » 3 4 *4 1/2. ) 

M2*M»EJ{ T* H/2.# 41#Pl/? ,* 31*01 /?. »x 2*P 1/2. # 12 « 3 1/2.# 4 3 
l*Tl/2.#3J« Ul/i.# A4 UU 2 . # 34 Ml/ 2 .) 

32*H*E4(T*H/2.,41*Pl/2.#3l«:U/2.#4?«!il/2,#P2#$l/2.#43 

1* C / !.# j 3* Ul/2 .* A4« Yl /2. *3 4# 4 1/?.} 
T2»"i*'5<r*H/2.*Al # *l/2.#31*3l/2»#42«"l/2.#3 ,< il/2 # *Ai 
l*Tl/2.*3 3*Lt/2.#AA*Xt/2.#34«4l/2.) 

U2*4*£iCT* H/2., Al4Pt/2.,3i«Ql /?.# A»45i/2 -# .j2« 
l* U/2.#3i« Ul/2. #44* VI n. *34<« 1/2. ) 

Y2*M*£7(T4M/?.*4t*Pt/2.#31«Qt/2*#42«<M/2.#-324>l/2,#A3 
1« 11/2. #3 3* Ul/2.*44*tfl //.#3 a«4 1/2. > 

4<*H*£8( r. H/,;.# 4l*pt/2,,3l«01/2..A2*»l/2.,3 2« 51/2 .#a 3 
l* U/2.# J Ul/2 .# 14* YU 2. #34*41/2.) 

C 

?J*3*cl<r.,-i/2.# 4t#P2/2..3l«D?/2,^A2*P2/2.*3;2«S2/2.>x3 

U 12/2.# J I* U 2/2. . (4* 4>t 2. *3 4* 4 2/2.) 
ii*H*£2( T* H/2.* 41 *P2/ 2. # 31 ♦ 42/2.* 4 248 2 / 2 ., 3 J* $ 2/». #* 5 





I« !?/>•» 33 * U 2 / 2 * » A 44 V 2 / 2 *# 34 *W 2 / 2 . > 

* JM»£ 3 tT*H/w.* 4 i*P 2 /?.# 3 l* 0 l/ 2 .*A 2 « 42 / 2**324 $ 2/?., 4 3 
l*T 2 / 2.#334 U 2 /»,#*fc#V 2 / 2 ,#S 4 *'< 2 / 2 ,> 

SlM*£4CT«H/Ji»M' P?/?.#3l«42/2*#42«42/2 .p*?*$2/>»*43 
i*! 2 / 2 .* ijt , 44 *¥ 2 / 2 ** 344 * 2 / 2 . > 
r 3 *H**- 5 (T*H/ i:«P 2 / 2 ., 3 l 432 / 2 .,* 2 «« 2 / 2 ., 32 « 52 / 2 # ,*3 

1 412 / 2 ** 3 J< U 2 />. ,a 4 *¥ 2 / 2 .* 344 W 2 / 2 .) 

IpJ"'**£ 6 <r« H/ 2 .» 41 *P 2 /?.* 3 1 * 72 / 2 .^ 2*4 2/2 , , 32 * 32 /’ ,# 4 3 
l« 12 / 2 **334 U 2 / 2 .* 44 4 * 2 ' 2**34 4 42 / 2 * ) 

41 *P 2 / 2 .» 31 * 32 / 24#4 2 *P 2 / 2.#5 2 ts 2 / 2.»4 3 
l»? 2 / 2 ** 33 *;,?/ 2 . ' 44 *V?/ 2 ,» 34 <* 2 / 2 .) 

«M*£ 3 {N H/ 2 • « «l <P 2 /J .#81 «!) 2 / 2 .» *»♦* 2/2 .# }?« $ 2/*.*4 J 
I* 12/ 2 .#3 ?♦ U 2 /?.* 44 « 42 / 2**34 44 2 / 2 * ) 

P 4 *H*£i(T»H* 4 t«pJ# 3 i«Q 3 * 4 2 *c U 3 2 * $ 1 #4 J* X 3 # J J* U 3 *4 4 * 43 
1 * 2 «** 4 3 ) 

34 «^*C :CT*^*A|*P 3 , 81 * 2 3 * 42*4 3*9 2*5 3*4 3 * T 3 * 3 3 # U S*A 4 *v *3 
l * £4* 4 J ) 

S<.*H*£ 3 (T«H. 4 l«PJ# 3 !* 33 # A 2 *r * ! , 9 2 « i J, 4 J* T 3 # J %*„ 3*4 *♦ *3 
t**i44J) 

* 4 «*«£ 4 C f* H» 4 |*»J» 3 l 4 )S ,4 2 " J * 32 * 3 3* 4 5 * T 3 * 3 J *'J 3 * 44 * 3 
i * £4 • 4 3 ) 

r 4 « H*£ 5 (T ♦ Hp 41 4 P 3 * 31 *^ 3 # *2 «P 3 . 32 < 33 * 43 * Ti ,3 34 U 3 * 44 * «# 3 

1 ** 4 * 43 ) 

w.<.sH*£4(T«H*4l4pJ**l*j3,4 24F$* ! 32*S3*4 34 TS*3J4„3 i ,44*YJ 
l » i 4 ♦ '4 i ) 

/ 4 a 4 *£ 7 f T« Hp 4 i«P 3 » 31 «‘) 3 » 4 2«4 . 3*3 2 « j 3*4 3 * 1 S »9 3 *d 3 fi 4 « *3 

1 * *4 * 4 3 ) 

44»:*#£J« T« ii.*i #P3»81 433 , 4 ?4fi i*3 2*5 ip 4 J* n#3 3« i JSf44«/3 
1 * j* ♦ 4j J 

A 1« 4 i 4 ( (»l *? ,*»2«2, *P J 404 )/§, ) 
Jl*U 4 <oi« 2 '*i 2 f 2 ** 1 J«;> 4 )/ 6 . , 

4<»4?4((^j J* 2 . *^ 3 * 1 ; V)/,, , 

• 12 «* 2 *<(> 14 ?.» > 24 *.*s j« 34 )/ 6 . ) 
4i»U4(CTl*?.*T2«2.*Tj4T4)/*. ) 
i J«M« ((Ui * 2 .*U! 42 . 4 UJ« U 4 >/ 3 . J 
» 4 *i 4 p((Vl«?.*rf 2 *?.*V?*V 4 )/b, ) 
3;*j4*((^i«?.*V242,*tJ**4)/4. > 

T»T« 4 


f , HJ)* 4 l *:3 3 (r) 4 ji 4 ' 5 lV{T) 

31 ( J >«A 2 *C 0 S{T )* J 2 * 5 H(T ) 
*’2CJ»»4i*C0S<>.#T)*334S!S<(2.* t) 

G2( J>*44»: <J3( 2. *1)4 j«*sj ’4( 2** n 

L*IM 

IF (L •£}* 3 ) JO TO 113 

*c ro 1 -jo 

lly 4FIT£ <4*3C > T p 41# i t » A 2# * >, 4 J ,J | * 4 4 ** 4 


50 FC'-'MI (1X#F2. 4# 3x*Fl 0.5#2X.F 10.5#3x#F 13 .5#2X #F 13 .5 #3X 
1>F13.5#2X»F13.5#3X#F10.5#2X#F10.5> 

L *1 

i - C CC^T I Nu£ 

L*1 
T = 0 

HMTE (5»U> 

41 F3MAT (IX# ///# IX* 3X# *TI • » >X# ‘Fi* #1 OX# • F2* > llX# 1 Gl* 
1# 12X# *G2*# / ) 

C 

OQ Z: C J = 1 * N M A X 
T=r« j 
L=L« 1 

IF (L .£0. 3) GO TO 210 
GC TO 200 

210 -I FIT: ( o * 5 0 0 ) T#F1(J)#F2(J)#G1(J)#32< J) 

500 FO*- < \ T Cl X#F2.4* 3X#F1 0.5#2X,? 10.5# 3X#F10. 5#2X»F 10.5) 

L = 1 

20 0 CC-iTIN'wi 
C 

CALL t^ir 
EM 


#-\t* t r*- TKl \ 1 1 


PAGSi & 


iwnwiinr 


APPENDIX D 

PROGRAM FOR EXACT SOLUTION OF 
STANDING WAVE CASE 
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C G»'*Y H MCDONALD 

C analytic solution-standing wave case 

DIMENSION F I ( 5 20 )# F 2( 5 OO) 

<» 1 = 3 . 14159 26 
H*D. 1 
X »=0 • 0 
MNAX=50Q 
E F 3~ 1 • l 
W9AR = C « l 
L = i 


?£AD IM GAS DYNAMIC INDEX AND INTERACTION INDEX 

F £ AD (5# 5) I#XN 
5 FORMAT (I5»F1C.5> 

INITIAL CONDI TIQNS 

5 E AD (5*1) ) C 1 » C 3# PHI 1 » PHI 3 
1C FORMAT (4F10.4) 

WRITE ( & » 1 5 ) I 

15 FCRMAT (lx# * GAS DYNAMIC I NOE X • » I 5) 

WFITE (S#1S > X M 

U FORMAT UX# ‘INTERACTION INDEX *»FF»2#///) 

WRITE C$#20) C1#C3#PHI1#PHI I 

23 FORMAT (IX#* INITIAL CONDITION S'# / /#1X#4X# 'Cl* »F ID • 4 # / 
i . iX»4X# , C3 , #Fl0.4#/#lX»4X»'PHIl , #Fl,Q.4»/»lX»4X»‘ o HI3* 
Z»fl0.4#/) 

WFITE C » » 3 C ) 

30 FORMAT ( IX »?X# * T IMF* , 14X# * Al* # 1 3 x# » A 3 ♦# l 3X# • 3 l* #12X 
1# *3 3* #15<» »CI‘#12X,‘C3 '#1 IX,* PHI l SIDX# 'PHI 3‘#/) 

JO IDG J = 1 # N M A X 
X-=X J4£°S 

IF (XN .ED. D..AND. I .EQ. 1) GO TO 6 
IF ( XN .ED. 40. .AND. I • E3. 0) GO TO 50 

\C GAS DYNAMICS 

50 wCNT I \ jE 
S « = I . J 

S=( 3 5353 *XN*C PS *(-!.)*• J) 

T*I. -:x°(-0.5*w8 AK* XJ) 

U - 3 * T 

/»!. /cost u > 

w=EXPC*0.5*w3AF* XJ) 

C1=W** 

Z = (S IN(U)/COSC U) ) 

C3 = ( W/2.J234 )-l 
PhIl=PHIl 
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phi 3 * 2 - **»h i i -< )/ 2 . >*p i 

sg n 200 


C NO CON 3 US TI ON 
C 

So CCNTIaLE 
*,=2. 

:S*(2P5*(-1.)**J>/C2.*X9AR> • 

CT=1 . “SXPC-C.5* w34R*XJ> 

Cl = C S *CT 

C3=2. /(EXP (CUM l,/SXP(CU>> 

C*=EXP (- 5. 5* m3 AR *x J ) 

Cl =C 

C2*CEXPCU )-l./ZXP(CJ> > /(£XP< CU>« 1. /• XP( cu> > 

:*=t:/«/ 2 .) *cz 

Phi 1=PHI1 

PHI3-2.* 3 HI I-? J*PI 
20C CONTINUE 

A 1 = C 1 * COS ( P M 1 1 ) 

31=C1*G INCPHII ) 

A 3 a G 3*C3S( pH ! 3 ) 

3 5 - C 3 * 5 I \' ( PHJ3 ) 

Fit J )=A1*CJSCXJ )«-31*S!X( XJ) 

F2( J ) = A 3*CQ SC2 . * XJ>«:3 3*SI*N( 2. *XJ ) 

IF CL .E3. 2) GO TO 5 0w 
L = L* I 
GO TO 1'JJ 
500 T I i 2 a X J 

*iF I T £ (6* 5000) T I M E* A I# A 3*3 1 *33*C1*C 3 »P,-| I 1 * PH I 3 
5 000 FCP‘U T C1X»F7.6, 3X * FI 0. 5* 5 X*F 10*5* 4 X, f 10. 5# 4X*F10.5*SX 
1 * FI .’.5*5X*F1 0.5, 4X*F 10.5, 4x*- 10.5) 

L = i 

100 CONTINUE 
C 

c 

x„ = o. 

L = t 

WrlTE (5*5000) 

6000 FQPN A T ( i X , ///, l X, 3 X* * Tl ME » * 1 2X * »F 1 » » 1 5< » »F 3 * ,/ ) 

C 

CC 300j = 1,NM4X 
X„ = .< J«H 

! F (L .£3. 2) GO TJ 700 
L~L* l 
GO TO 300 
700 TH2 = XJ 

*RIT2 (6*7000) TIM£*F1(J),F2( J) 

7000 FORMAT C 1 X ,F 7 . 4 * <3 X *F 1 0 .5 » 7 X »• 1 0. 5* / ) 

L= 1 

30 C CONTINUE 


£ 


^Z u l j a ^ o o ITRE 


poor 
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C 3 AF Y H MCOQNALQ 

C ANALYTIC S OLl* T 1 ON ■ T F A Y £L L 1 N 3 *AY£ CASE 

OH* M SION fl(50 0)»F2(503)*jl( 500 >»3 2( 500) 

?I 8 J. 1415926 
H s 0 • 1 
X»s«).0 
NM4X 8 5O0 
2 F S 8 0 • 1 
W34* 8 0.l 
L*1 
C 

C P t AO IX 3 AS DYNAMIC IN02X AND IN T£* ACT ION INDEX 
C 1 

S£ AO (b#5) I.XN 
5 F3FMAT (X5»F10.5) 

C 

C INITIAL CONDITIONS 
C 

n£AD (5*13) Al»A2sA3»A4 
1C FCAMAT (4F10.H ) 

*hIT£ (6»15) I 

15 rCFMAT (IX# •GAS DYNAMIC 1M0£X*»I5) 

W c I T - (6 >1 6 ) XN 

to FOAM A T (lX,M.4feSACTI3N I .\IQ£X * >F 7 • l '/ // ) 
mAH l (5,20) A 1 , A2»A3/A4 

20 FCMAT UX»*1n1T1AL CON 0 1 TIG'! S • * // * IX * 4 X , *4 1 * ,f 10 . 4 * t 
l#lX»AX#*42 , »Fn*4»/#iX*4X#*4j*#FlC.W*lX*4X»*A4i 
2# FI 0. 4»/) 

WAIT- (S»iO) 

10 FCrfiAT ( IX#4X# *TIM£*» 14X* *41* *15X»' A2»* 15X#*43* »17X 
1# *44 »#/) 

CO DO J 8 l*.NMAX 

X^-XJ4£pj 

IF ( XN .El. 0..4YD. ! .£3. 1) GO TO 6 0 
IF CXN .EG . 2).. AND. I .E3. 0) GO TO 5 3 
C 

C NO G A 5 DYNAMICS 
C 

50 CCMTINLE 

5 *0 . 7 07 11 * £?S* XN 
T*l. -EX»{-D.5*W3Ar *XJ> 

0 = S*T 

•/-l. /CO SC J) 

MsEX=»(-'j.5*W3A-i*XJ) 

42 = ** V 

i*3lN(U>/C0SCJ) 

A 4 = ( W / 2 . 6 .2 6 4 Z ) * 2 
GO TO 20 j 
C 

C NC CC^SUSTION 
C 


■*4 


'A - 


St) CONTINUE 
C$*£PS/WJAP 

CT*1«“£XP(“w«5*h3AR*XJ) 

CU*C$*CT 

Ct/»2. / (EXP (CU)« l./EXPCC'J>> 

C* = EX?(-0. 5*w3AR*XJ> 

A 1 * C X * C / 

C2*(EXP(CU)-1./EXP(CU) >/(£XP( CUXl. /EXP(Cu)) 
A!=(C*/2.)*CZ 
3C TO 200 
200 31»-A2 
3 2 = A 1 
33= -A A 
3A* A 3 

Fl( J)*A1*:0S(XJ)*31*S1N( XJ) 

G1CJ )*A2*C0S(XJ )432*SIN(XJ) 

F2< J )= A 3*0 OS (2. *XJ)«3 3*5I.N(2. *XJ> 

32( J>*AA*C0SC2.*XJ)4-3A*SIN(2.*XJ) 

IF (L .E3. 2) 00 TO 303 
L*L* l 
GO TO 103 
500 THE = XJ 

<tf>I T E (5* 50 03) TIME*A1» A2*A3*AA 
5000 FCMAT UX*F7.A*7X* (A(F12. 6*5 X)>) 

L = l 

100 CON T In U t 


X„ = 0. C 
L *1 

WFITE (5*5033) 

6000 FO-'MI (1 X* ///* IX* 3X* TI M£ • *1 2X* *Fl •» I5X » »F2* *16X» '31 ' 
l* IAX* *32'* / ) 

00 30 C J = 1 » N 4 A X 
XJ=XJ4 H 

IF (L .£3. 2) 30 TO 700 
L = L* l 
3C to 300 
700 TIME*XJ 

AH1TZ (6*7000) TIME*FUJ)»F2( J)*0l( J)*G2( J) 

7 OtO F QM A T < IX *r 7 . 4 , 3X * F 1 0 . 5 * 7 X *" 10 . '5 *6 X * F 10 .5 » 7 X *F 1) . 5 ) 
L*i 

30C CONTINUE 
CALL EXIT 
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PRESENTATION OF ACOUSTIC 
PRESSURE CALCULATIONS 
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ACOUSTIC PRESSURE DERIVATION 

To calculate expressions for acoustic pressure, recall equation 
(2.48) which stated 

-(« ||- + h e 2 (u 2 + 2u ^ 

P - P * e \ V 3z // . (F.l) 

This equation represents the unsteady state deviations of acoustic 
pressure. When expanding equation (F.l) into a Taylor series expansion, 
the resulting equation becomes 

p = p = i . E |i + e 2 [jj(G 2 + h • h) -5 1| + h (H-) 2 ] 


+ . . . . (F.2) 

Recall that the steady state solution was represented in equation (2.35) by 



When expanding (F.3) into its Taylor series expansion, the result becomes 


P . 1 - *se 2 fei) 2 
\ dZ / 


+ • . 


(F.4) 
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where p is the steady state acoustic pressure. Therefore, the difference 
in general acoustic pressure and steady state pressure can be expressed 
by subtracting equation (F.4) from (F.2). For this investigation, a 
restriction on the velocity potential 4 was that it was a function of 6 
and t only. In doing this, the pressure difference equation becomes 


P - P 


-e 



e* 

T 



(F.5) 


Using the same Fourier series expansion for the velocity potential 4 as 
expressed in equation (3.18), the acoustic pressure difference equation 
(F.5) can be expressed in terms of the product of modal amplitudes and 
trignometric function in the transverse 0 direction. Substituting the 
appropriate forms of equation (3.18) into equation (F.5) and simplyfying, 
the resulting pressure difference equation become 
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Since the coefficients In equation (F. 6) are functions of tine only, 
these coefficients have been included in the calculations of the program 
in Appendix B. Thus, for i.;iy given angle 8, values for the modal 
amplitude at any given time range can be calculated therefore determining 
the acoustic pressure difference of that desired location. 


